CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > General Forums > Main CFD Forum

Numerical treatment of the source term in combustion equations

Register Blogs Community New Posts Updated Threads Search

Like Tree10Likes

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   September 7, 2016, 16:37
Default
  #21
Super Moderator
 
Tobi's Avatar
 
Tobias Holzmann
Join Date: Oct 2010
Location: Tussenhausen
Posts: 2,708
Blog Entries: 6
Rep Power: 51
Tobi has a spectacular aura aboutTobi has a spectacular aura aboutTobi has a spectacular aura about
Send a message via ICQ to Tobi Send a message via Skype™ to Tobi
Hi all,

well I want to split the calculation into:

a) Calculate the source terms omega and update the Y values
b) Solve the transport equation

At the moment I am only interested in the steady-state solution (no transient).

for a) I would do the following:

[LaTeX Error: Syntax error]

Solving that stiff system with the BDF like:
\frac{Y_k^n - Y_k^{n-1}}{\Delta t}= S_k^n

Y_k^n = S_k^n \Delta t + Y_k^{n-1}

Here I update the "old" values with the BDF method and then I solve the transport equation (simple one). Of course the above equation has to be looped till I get steady-state here.

Quote:
Typically, for non linear coupled equations you still need the jacobian to determine the maximum allowable time step. For each equation you get the maximum absolute eigenvalue of the jacobian and the max dt is something inversely proportional to that.
Ahhh okay. So after I have the Jacobian and the Eigenvalues (Eigenvalue calculation is already implemented), I can estimate the time-step for the chemical reaction like the inverse of the largest Eigenvalue? Good to know.
__________________
Keep foaming,
Tobias Holzmann
Tobi is offline   Reply With Quote

Old   September 7, 2016, 16:47
Default
  #22
Senior Member
 
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 6,768
Rep Power: 71
FMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura about
Actually, if you consider the splitting and you adopt an implicit method, you have to consider that, say Y* the intermediate field, the equation writes as


Y*- dt S(Y*)=Y^n-1
FMDenaro is offline   Reply With Quote

Old   September 7, 2016, 19:43
Default
  #23
Senior Member
 
Michael Prinkey
Join Date: Mar 2009
Location: Pittsburgh PA
Posts: 363
Rep Power: 25
mprinkey will become famous soon enough
Quote:
Originally Posted by Tobi View Post
@ Raj ... wrong thread.



I meant the Backward Differentiation Formula of step size 2 (BDF 2) - like here: https://en.wikipedia.org/wiki/Backwa...iation_formula

At the moment I am implementing the BDF1 that corresponds to the Backward Euler discretization.
I think that you want Adams Bashforth Moulton, rather than backward differencing for higher-order time integration.

https://en.wikipedia.org/wiki/Linear_multistep_method

ABM schemes use the RHS value at two or more different time levels, rather than using multiple old time values of the solution. They generally work better than high-order backward differencing. If you are a building a split implicit/explict scheme, Trapzoidal Rule (2nd-order Adams Moulton) for the implicit terms and 2nd-order Adams Bashforth for explicit terms often work well.
FMDenaro and Tobi like this.
mprinkey is offline   Reply With Quote

Old   September 8, 2016, 03:13
Default
  #24
Super Moderator
 
Tobi's Avatar
 
Tobias Holzmann
Join Date: Oct 2010
Location: Tussenhausen
Posts: 2,708
Blog Entries: 6
Rep Power: 51
Tobi has a spectacular aura aboutTobi has a spectacular aura aboutTobi has a spectacular aura about
Send a message via ICQ to Tobi Send a message via Skype™ to Tobi
Quote:
Originally Posted by FMDenaro View Post
Actually, if you consider the splitting and you adopt an implicit method, you have to consider that, say Y* the intermediate field, the equation writes as


Y^*- \Delta t S(Y^*)=Y^{n-1}
Yes you are right, I should write some ^* to indicate a intermediate state. Just to clear my mind,... if we say splitting we mean, that we solve the source term first and then the transport equation.

Quote:
I think that you want Adams Bashforth Moulton, rather than backward differencing for higher-order time integration.
I will check it out. Thanks for the information.
__________________
Keep foaming,
Tobias Holzmann
Tobi is offline   Reply With Quote

Old   September 8, 2016, 04:16
Default
  #25
Senior Member
 
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 6,768
Rep Power: 71
FMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura about
Actually, you cannot get higher order accuracy even using AB methods, the reason being in the splitting. Immagine you work first without splitting. Integrating in time you get

Y^n+1-Y^n = Int[tn+1;tn] ( D d2 Y/dx^2 + S(Y)) dt

Now, with the splitting you solve

Y*-Y^n = Int[tn+1;tn] S(Y*) dt

Y**-Y* = Int[[tn+1;tn] ( D d2 Y**/dx^2) dt

Therefore, the best you can get (but with some care), is only second order accuracy in time.
FMDenaro is offline   Reply With Quote

Old   September 8, 2016, 06:16
Default
  #26
Senior Member
 
sbaffini's Avatar
 
Paolo Lampitella
Join Date: Mar 2009
Location: Italy
Posts: 2,151
Blog Entries: 29
Rep Power: 39
sbaffini will become famous soon enoughsbaffini will become famous soon enough
Send a message via Skype™ to sbaffini
Still, unless there is some specificity regarding chemical source terms, as he intend to just reach the steady state, i guess such effort on the time acuracy is not actually worth it. He jest needs to keep it stable and consistent, not time accurate.
sbaffini is offline   Reply With Quote

Old   September 8, 2016, 06:28
Default
  #27
Senior Member
 
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 6,768
Rep Power: 71
FMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura about
That's true Paolo, however at the steady state only the time derivatives in the truncation error vanish. I am not sure If the remaining spatial derivatives retain some time step multiplication so that the temporal accuracy enters into the global accuracy.
But I agree, I wouldn' take too care of time accuracy in a first stage of a steady solution.
FMDenaro is offline   Reply With Quote

Old   September 8, 2016, 06:34
Default
  #28
Senior Member
 
sbaffini's Avatar
 
Paolo Lampitella
Join Date: Mar 2009
Location: Italy
Posts: 2,151
Blog Entries: 29
Rep Power: 39
sbaffini will become famous soon enoughsbaffini will become famous soon enough
Send a message via Skype™ to sbaffini
Quote:
Originally Posted by FMDenaro View Post
That's true Paolo, however at the steady state only the time derivatives in the truncation error vanish. I am not sure If the remaining spatial derivatives retain some time step multiplication so that the temporal accuracy enters into the global accuracy.
But I agree, I wouldn' take too care of time accuracy in a first stage of a steady solution.
That should not happen with the method of lines and the same should apply also for the time splitting which, if i am not wrong, is not practically different from a segregated approach.
sbaffini is offline   Reply With Quote

Old   September 8, 2016, 06:45
Default
  #29
Super Moderator
 
Tobi's Avatar
 
Tobias Holzmann
Join Date: Oct 2010
Location: Tussenhausen
Posts: 2,708
Blog Entries: 6
Rep Power: 51
Tobi has a spectacular aura aboutTobi has a spectacular aura aboutTobi has a spectacular aura about
Send a message via ICQ to Tobi Send a message via Skype™ to Tobi
Dear all,

Quote:
Originally Posted by sbaffini View Post
Still, unless there is some specificity regarding chemical source terms, as he intend to just reach the steady state, i guess such effort on the time acuracy is not actually worth it. He jest needs to keep it stable and consistent, not time accurate.
Thats true. At the beginning I will not investigate to deep into the time accuracy because I am only interested in the steady state situation. After that, I may implement a better time scheme but up to know I should first get a steady-state solution. Sometimes I have some mind blockers if I start implementing the calculation (which field has to be updated and are not allowed etc). In any case the discussion with you is very helpful and interesting.

Quote:
I think that you want Adams Bashforth Moulton, rather than backward differencing for higher-order time integration.

https://en.wikipedia.org/wiki/Linear_multistep_method

ABM schemes use the RHS value at two or more different time levels, rather than using multiple old time values of the solution. They generally work better than high-order backward differencing. If you are a building a split implicit/explict scheme, Trapzoidal Rule (2nd-order Adams Moulton) for the implicit terms and 2nd-order Adams Bashforth for explicit terms often work well.
Very interesting, the last days I learned so much about the time integration ... I think it is not a "too big deal" to implement such a integration method but up to know I keep the BDF (1) for conservation and simplicity.

If you are interested, I will keep you up to date with the latest results. Unfortunately I will not find time the next 4 days ...
sbaffini likes this.
__________________
Keep foaming,
Tobias Holzmann
Tobi is offline   Reply With Quote

Old   September 12, 2016, 16:10
Default
  #30
Super Moderator
 
Tobi's Avatar
 
Tobias Holzmann
Join Date: Oct 2010
Location: Tussenhausen
Posts: 2,708
Blog Entries: 6
Rep Power: 51
Tobi has a spectacular aura aboutTobi has a spectacular aura aboutTobi has a spectacular aura about
Send a message via ICQ to Tobi Send a message via Skype™ to Tobi
Dear all,

it would be very helpful to get some feedback about section 1.2.2 in the attachment. The PDF is just a draft so nothing is finished

http://holzmann-cfd.de/forumCases/UserGuide.pdf
__________________
Keep foaming,
Tobias Holzmann

Last edited by Tobi; September 12, 2016 at 16:32. Reason: Corrected section 1.5 to 1.2.2
Tobi is offline   Reply With Quote

Old   September 12, 2016, 16:19
Default
  #31
Senior Member
 
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 6,768
Rep Power: 71
FMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura about
Quote:
Originally Posted by Tobi View Post
Dear all,

it would be very helpful to get some feedback about section 1.5 in the attachment. The PDF is just a draft so nothing is finished

http://holzmann-cfd.de/forumCases/UserGuide.pdf
there is no section 1.5
FMDenaro is offline   Reply With Quote

Old   September 12, 2016, 16:29
Default
  #32
Super Moderator
 
Tobi's Avatar
 
Tobias Holzmann
Join Date: Oct 2010
Location: Tussenhausen
Posts: 2,708
Blog Entries: 6
Rep Power: 51
Tobi has a spectacular aura aboutTobi has a spectacular aura aboutTobi has a spectacular aura about
Send a message via ICQ to Tobi Send a message via Skype™ to Tobi
Oh I am sorry, ... I ment 1.2.2.
__________________
Keep foaming,
Tobias Holzmann
Tobi is offline   Reply With Quote

Old   September 12, 2016, 21:41
Default
  #33
Senior Member
 
Arjun
Join Date: Mar 2009
Location: Nurenberg, Germany
Posts: 1,272
Rep Power: 34
arjun will become famous soon enougharjun will become famous soon enough
Tobias

What is the size of system we are talking here. Does it depend on number of species we are talking about.

So if you took the coupled system then does the system is of size Nspecies x Nspecies.

If thats the case can't you use direct solver on coupled system, you can use some multifrontal solver. I believe some open libraries are there that provide it.
arjun is offline   Reply With Quote

Old   September 12, 2016, 22:59
Default
  #34
Senior Member
 
Michael Prinkey
Join Date: Mar 2009
Location: Pittsburgh PA
Posts: 363
Rep Power: 25
mprinkey will become famous soon enough
I don't completely understand your flamelet formulation, but I grasp the essential ideas you are working towards. I am a bit concerned that you are choosing a timestep for the X updates and then iterating to convergence while picking a different timestep for your Y updates. This *may* be OK in the limit of large time, but my gut tells me that it does not preserve the relative time scales of the reaction and the transport terms. This convergence loop seems like it should essentially just give you equilibrium chemistry results at every timestep. I think I can follow your recipe and choose the two different timesteps independently and generate very different physical results and, of course, that indicates that the scheme in inconsistent.. Again, I haven't proved or tested this, but I see it as a potential problem.

Stiff (split) chemistry solvers usually assume a transport time step (DT) based on CFL or some other criteria. Then, reaction timesteps (dt) are chosen to evolve the system from t_0 to (t_0 + DT). Usually, dt << DT. But, you step forward the reactions to t_0 + DT, and then you evolve the transport fractional step to t_0 + DT. This retains consistency between the transport and the reaction timescales automatically and gives time accurate results. This will, of course, provide consistent steady state results as well as you take t->inf. But choosing different timesteps changes the effective timescale of each--techniques exist to precondition pseudo-timestepping schemes that discard time accuracy to get faster convergence, but those techniques are complicated and maintain the relative scale of the terms as time goes to infinity.
mprinkey is offline   Reply With Quote

Old   September 13, 2016, 02:04
Default
  #35
Super Moderator
 
Tobi's Avatar
 
Tobias Holzmann
Join Date: Oct 2010
Location: Tussenhausen
Posts: 2,708
Blog Entries: 6
Rep Power: 51
Tobi has a spectacular aura aboutTobi has a spectacular aura aboutTobi has a spectacular aura about
Send a message via ICQ to Tobi Send a message via Skype™ to Tobi
Quote:
Originally Posted by arjun View Post
Tobias

What is the size of system we are talking here. Does it depend on number of species we are talking about.

So if you took the coupled system then does the system is of size Nspecies x Nspecies.

If thats the case can't you use direct solver on coupled system, you can use some multifrontal solver. I believe some open libraries are there that provide it.
Hey,

the number of systems depend on two things. First the space discretization (mixture fraction) depend on the number of points we want to have. Lets say 100 and within one discrete point, the chemical source term depend on the number of species that are included within the detailed chemistry. The Jacobian is a nSpecies x nSpecies matrix.

The fact that I do not want to use some open libraries are as follow: I want to make my hands dirty to understand things in more details I know that it is somehow stupid to re-invent the tire but for a good understanding, it is necessary for me (I am not a genius - unfortunately)

Quote:
Originally Posted by mprinkey
I don't completely understand your flamelet formulation, but I grasp the essential ideas you are working towards. I am a bit concerned that you are choosing a timestep for the X updates and then iterating to convergence while picking a different timestep for your Y updates. This *may* be OK in the limit of large time, but my gut tells me that it does not preserve the relative time scales of the reaction and the transport terms. This convergence loop seems like it should essentially just give you equilibrium chemistry results at every timestep. I think I can follow your recipe and choose the two different timesteps independently and generate very different physical results and, of course, that indicates that the scheme in inconsistent.. Again, I haven't proved or tested this, but I see it as a potential problem.
I understand, it was the first thought that I had. Like:

  • Calculate a time step for the transport \Delta t_tr based on the Fourier number (I have no convection and therefore no CFL, hence, I calculate the time step based on the Fo number; for 1D calculation Fo < 0.5)
  • Go into the calculation for the source term; here, choose another timestep \Delta t_ch and go on in time till we reach \Delta t_tr -> consistency in time
Quote:
Originally Posted by mprinkey
Stiff (split) chemistry solvers usually assume a transport time step (DT) based on CFL or some other criteria. Then, reaction timesteps (dt) are chosen to evolve the system from t_0 to (t_0 + DT). Usually, dt << DT.
I will change the schematic (which is still not good) today evening.

Quote:
Originally Posted by mprinkey
But, you step forward the reactions to t_0 + DT, and then you evolve the transport fractional step to t_0 + DT. This retains consistency between the transport and the reaction timescales automatically and gives time accurate results. This will, of course, provide consistent steady state results as well as you take t->inf. But choosing different timesteps changes the effective timescale of each--techniques exist to precondition pseudo-timestepping schemes that discard time accuracy to get faster convergence, but those techniques are complicated and maintain the relative scale of the terms as time goes to infinity.
Okay, then I will change the stuff to a consistency behavior. Seems to me a better choice too.

Note that till now I am not solving the temperature equation. I hope to find time to add this to the figure in section 1.2.2 today evening.

Thanks for your replays.
__________________
Keep foaming,
Tobias Holzmann
Tobi is offline   Reply With Quote

Old   September 13, 2016, 15:08
Default
  #36
Super Moderator
 
Tobi's Avatar
 
Tobias Holzmann
Join Date: Oct 2010
Location: Tussenhausen
Posts: 2,708
Blog Entries: 6
Rep Power: 51
Tobi has a spectacular aura aboutTobi has a spectacular aura aboutTobi has a spectacular aura about
Send a message via ICQ to Tobi Send a message via Skype™ to Tobi
I made changes to the calculation procedure. I really appreciate your opinions. Thanks in advance. The file is updated and is the same as before, or just »click«
__________________
Keep foaming,
Tobias Holzmann
Tobi is offline   Reply With Quote

Old   March 20, 2017, 08:58
Default
  #37
Senior Member
 
Arjun
Join Date: Mar 2009
Location: Nurenberg, Germany
Posts: 1,272
Rep Power: 34
arjun will become famous soon enougharjun will become famous soon enough
Quote:
Originally Posted by Tobi View Post
Offtopic:

Hi Arjun, I made my own code that loads chemistry, transport and thermodynamic data. The code reads a dictionary where the user can specify the fuel and oxidizer stream and other solver related stuff. After that, I construct the polynomials for fitting (viscosity, thermal conductivity, binary diffusion etc.). However, the user can specify the calculation of the gas kinetic equation and if he want to use fitting procedure or direct calculation. Then I build my flamelet class with all information which then should solve everything. Everything is already implemented, TROE, LOW, ENHANCED, SRI ... but the only thing that is missing is the solving of the chemistry. I was thinking of using sub-time stepping, like flamelet equation dt and then solve chemistry with an ODE in order to reach that dt. After that go to the next time step (I think it is a common approach). At the moment (or 3 weeks ago) I started to implement simple Euler Integration which is easy but not really powerful. However, doing that, I would have the basis of other solver implementations. But like always - no time. So the project is stlil at the moment. If you are interested: http://afc.holzmann-cfd.de/ but it is outdated. If you want, I can upload the latest doxygen stuff today evening.
Okay I made progress in my solver. Now i am able to run full species transport equation along with reacting multicomponent fluid with constant diffusivity. (I need to check why solver acting up with diffusivities varying, so it is next thing to do).

The current version uses VODE based species stiff ode integration and I use operator spliting method that fluent and starccm typically offer.
StarCCM also offers something called Chemistry ADI (This is what i believe they call it) that can handle even larger time steps. I am not looking at it at the moment but in future if i am unhappy with operator splitting i might visit it.

Next week or so I try to validate the solver with various sets up and see how robust it is. So far looks good but it has not done any challenging thing too.

Last edited by arjun; March 20, 2017 at 10:47.
arjun is offline   Reply With Quote

Old   September 15, 2020, 13:42
Default
  #38
Senior Member
 
Reviewer #2
Join Date: Jul 2015
Location: Knoxville, TN
Posts: 141
Rep Power: 10
randolph is on a distinguished road
Hi All,

I am creating a simple solver that involve transport and reaction. Different from Tobi's purpose, I would like to get the time accurate solution. I think this post is relevant and also extended from the original discussion.

Here is what I have done so far. I would like to hear your valuable advices to improve the solver's efficiency since I am really interested in simulate the transport and reaction for a large timescale ~ hr to day.


1. First I solve the reaction term dY/dt = f(Y) only with the ODE solver in OpenFOAM from t0 to t0+Dt with Dt being the time step size of the transport equation. Y0 and YEnd is the initial and final results at t0 and t0+Dt in this ODE.

The ODE time step size (dt) that is automatically determined by the ODE solver of OpenFOAM so that sufficiently small to meet the error tolerance.

2. Second I construct the source term as dYdt = (YEnd-Y0)/Dt.

3. Third shuffle the dYdt as an explicit source term to the transport equation and solve the transport equation.


Here are my questions:

I would like to take the time step (Dt) as large as possible. But I observe that the transient solution has depends on the time step size in the transport solver (Dt). The solution changes even I reduce the time step size to 1e-2, which is very challenge to the intended simulation length (hr to day).

What is the bottleneck of the code?
I guess this way to split the equation is only first order accurate. Will Strang splitting help? Without too much of work, Strang splitting can be implemented with the first order time integration scheme for the transport equation. Based on the exsiting first order splitting, I can use second order scheme for the transport equation. If one choose the accuracy between splitting and transport, which one' accuracy should be prioritized?

How can I estimate the time scale of the reaction?

Thanks in advance,
Rdf

Last edited by randolph; September 15, 2020 at 15:35.
randolph is offline   Reply With Quote

Reply


Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
what is swap4foam ?? AB08 OpenFOAM 28 February 2, 2016 01:22
[Other] Adding solvers from DensityBasedTurbo to foam-extend 3.0 Seroga OpenFOAM Community Contributions 9 June 12, 2015 17:18
OpenFOAM without MPI kokizzu OpenFOAM Installation 4 May 26, 2014 09:17
[swak4Foam] Swak4FOAM 0.2.3 / OF2.2.x installation error FerdiFuchs OpenFOAM Community Contributions 27 April 16, 2014 15:14
DecomposePar links against liblamso0 with OpenMPI jens_klostermann OpenFOAM Bugs 11 June 28, 2007 17:51


All times are GMT -4. The time now is 14:51.