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

Why decreasing time step makes diffusion more significant for diffusion equation?

Register Blogs Community New Posts Updated Threads Search

Like Tree1Likes

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   February 15, 2018, 14:29
Question Why decreasing time step makes diffusion more significant for diffusion equation?
  #1
Senior Member
 
Join Date: Oct 2017
Location: United States
Posts: 233
Blog Entries: 1
Rep Power: 10
TurbJet is on a distinguished road
I am trying to solve a simple diffusion equation

\partial_{t}u + a\partial_{x}u = b\partial_{xx}u

where a & b are some given constant, for the sake of simplicity, let's say they are both 1, and the spatial domain is [-1, 1] with periodic BC, the initial condition is numerical delta function u^{0} = \delta(x-0) = \frac{1}{\Delta x}.

Solving this by simply using a naive schemes:

u^{n+1}_{i} = u^{n}_{i} - a\Delta t\frac{u^{n}_{i+1} - u^{n}_{i-1}}{2\Delta x} + b\Delta t\frac{u^{n}_{i+1} - 2u^{n}_{i} + u^{n}_{i-1}}{\Delta x^{2}}

When I compare results from different time step sizes, the smaller the time step size, the more diffusive it appears to be: the wave, or the delta function, spreads out much faster in space (smoothed out) when advance in time, with smaller time step size. But according to the scheme, lowering time step size \Delta t will obviously reducing the 2nd-order derivative, i.e., the diffusion term, as well, thus it should have lead to less diffusion for smaller \Delta t.

So where does this conflict come from?

Thx!

Last edited by TurbJet; February 15, 2018 at 15:48.
TurbJet is offline   Reply With Quote

Old   February 15, 2018, 15:18
Default
  #2
Senior Member
 
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 6,896
Rep Power: 73
FMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura about
Well, there are some issues to discuss. First, look at the local truncation error of your discretization:

LTE=-0.5*dt*(a^2*d2u/dx^2-2*a*b*d3u/dx^3+b^2*d4u/dx^4) +...

and, I am right, you have a higher order term in the form dx^3/dt (try to develop the LTE). Therefore, for fixed dx and vanishing dt you should get an increasing in the error. This is due to the fact that you need convergence for dx and dt going both to zero.

Second, you cannot use a Dirac initial function to test such PDE as it is a non regular function.
FMDenaro is offline   Reply With Quote

Old   February 15, 2018, 15:30
Default
  #3
Senior Member
 
Join Date: Oct 2017
Location: United States
Posts: 233
Blog Entries: 1
Rep Power: 10
TurbJet is on a distinguished road
Quote:
Originally Posted by FMDenaro View Post
Well, there are some issues to discuss. First, look at the local truncation error of your discretization:

LTE=-0.5*dt*(a^2*d2u/dx^2-2*a*b*d3u/dx^3+b^2*d4u/dx^4) +...

and, I am right, you have a higher order term in the form dx^3/dt (try to develop the LTE). Therefore, for fixed dx and vanishing dt you should get an increasing in the error. This is due to the fact that you need convergence for dx and dt going both to zero.

Second, you cannot use a Dirac initial function to test such PDE as it is a non regular function.
I see your point. but I am not sure about your formula for LTE: it seems all terms are with \frac{\Delta t}{\Delta x}, which still leads to decreasing with smaller \Delta t. Maybe some minor mistake in it?
TurbJet is offline   Reply With Quote

Old   February 15, 2018, 15:42
Default
  #4
Senior Member
 
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 6,896
Rep Power: 73
FMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura about
Quote:
Originally Posted by TurbJet View Post
I see your point. but I am not sure about your formula for LTE: it seems all terms are with \frac{\Delta t}{\Delta x}, which still leads to decreasing with smaller \Delta t. Maybe some minor mistake in it? And could you tell me how you derive this formula?

I wrote only the term of O(dt), of course you have also terms of O(dx^2). The LTE is determined by substituting the Taylor expansion in time and space in your FTCS scheme. You can do it by hand, is a bit long to manage, or you can use Maple.

However, you cannot use the Dirac as initial condition and the study of the convergence requires to use dt/dx->0
FMDenaro is offline   Reply With Quote

Old   February 15, 2018, 16:04
Default
  #5
Senior Member
 
Join Date: Oct 2017
Location: United States
Posts: 233
Blog Entries: 1
Rep Power: 10
TurbJet is on a distinguished road
Quote:
Originally Posted by FMDenaro View Post
I wrote only the term of O(dt), of course you have also terms of O(dx^2). The LTE is determined by substituting the Taylor expansion in time and space in your FTCS scheme. You can do it by hand, is a bit long to manage, or you can use Maple.

However, you cannot use the Dirac as initial condition and the study of the convergence requires to use dt/dx->0
Taylor series! How stupid I am !

Anyway, I wrote out the formula for LTE only with leading order terms (probably some mistakes in those coefficients):

LTE \approx -\Delta tu_{tt} - \frac{\Delta x^{2}}{3}u_{xxx} + \frac{\Delta x^{2}}{12}u_{xxxx} = \Delta t(-u_{tt} - \frac{\Delta x^{2}}{3\Delta t}u_{xxx} + \frac{\Delta x^{2}}{12\Delta t}u_{xxxx})

so by lowering the time step size, the LTE for time derivative would not increase, but the LTE for spatial derivative will. So which means it will increase the spatial error, and it appears to be like diffusion.

Am I correct?
TurbJet is offline   Reply With Quote

Old   February 15, 2018, 16:19
Default
  #6
Senior Member
 
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 6,896
Rep Power: 73
FMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura about
As you see, if dx is taken constant while dt is diminished, you have both a dispersion and dissipation effects. However, you should still expand d2u/dt^2 using the PDE equation
FMDenaro is offline   Reply With Quote

Old   February 15, 2018, 16:24
Default
  #7
Senior Member
 
Join Date: Oct 2017
Location: United States
Posts: 233
Blog Entries: 1
Rep Power: 10
TurbJet is on a distinguished road
Quote:
Originally Posted by FMDenaro View Post
As you see, if dx is taken constant while dt is diminished, you have both a dispersion and dissipation effects. However, you should still expand d2u/dt^2 using the PDE equation
Yes, now everything is clear. But what do you mean by "expand u_{tt} using the PDE equation"? Do you actually mean by u_{xx} instead u_{tt}?
TurbJet is offline   Reply With Quote

Old   February 15, 2018, 16:29
Default
  #8
Senior Member
 
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 6,896
Rep Power: 73
FMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura about
Quote:
Originally Posted by TurbJet View Post
Yes, now everything is clear. But what do you mean by "expand u_{tt} using the PDE equation"? Do you actually mean by u_{xx} instead u_{tt}?
I mean to develop

d2u/dt^2= d/dt ( -a*du/dx+b*d2u/dx^2) = d/dx (-a*du/dt+b*d/dx(du/dt)) = ...
FMDenaro is offline   Reply With Quote

Old   February 15, 2018, 16:34
Default
  #9
Senior Member
 
Join Date: Oct 2017
Location: United States
Posts: 233
Blog Entries: 1
Rep Power: 10
TurbJet is on a distinguished road
Quote:
Originally Posted by FMDenaro View Post
I mean to develop

d2u/dt^2= d/dt ( -a*du/dx+b*d2u/dx^2) = d/dx (-a*du/dt+b*d/dx(du/dt)) = ...
Oh~~, and then plug it back to the LTE formula derived above to replace the u_{tt} term, right?
TurbJet is offline   Reply With Quote

Old   February 15, 2018, 16:38
Default
  #10
Senior Member
 
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 6,896
Rep Power: 73
FMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura about
Quote:
Originally Posted by TurbJet View Post
Oh~~, and then plug it back to the LTE formula derived above to replace the u_{tt} term, right?

yes, there is some arrangement to manage...
FMDenaro is offline   Reply With Quote

Old   February 15, 2018, 16:39
Default
  #11
Senior Member
 
Join Date: Oct 2017
Location: United States
Posts: 233
Blog Entries: 1
Rep Power: 10
TurbJet is on a distinguished road
Quote:
Originally Posted by FMDenaro View Post
yes, there is some arrangement to manage...
Gotcha. Thx very much. You have been really helpful !!
TurbJet is offline   Reply With Quote

Old   February 15, 2018, 16:42
Default
  #12
Senior Member
 
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 6,896
Rep Power: 73
FMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura about
Again, you see that using the Taylor expansion requires the solution to be regular...
FMDenaro is offline   Reply With Quote

Old   February 15, 2018, 16:49
Default
  #13
Senior Member
 
Join Date: Oct 2017
Location: United States
Posts: 233
Blog Entries: 1
Rep Power: 10
TurbJet is on a distinguished road
Quote:
Originally Posted by FMDenaro View Post
Again, you see that using the Taylor expansion requires the solution to be regular...
Um, theoretically speaking, yes. But when I actually ran it with discretized form of Dirac function, in this case, 2/dx, it did not show some wired behavior aside from strong diffusion....
TurbJet is offline   Reply With Quote

Old   February 15, 2018, 17:05
Default
  #14
Senior Member
 
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 6,896
Rep Power: 73
FMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura about
Quote:
Originally Posted by TurbJet View Post
Um, theoretically speaking, yes. But when I actually ran it with discretized form of Dirac function, in this case, 2/dx, it did not show some wired behavior aside from strong diffusion....
That's not exact. Your initial solution has still a spike in a point so that you cannot define the derivatives across it.
Also from a numerical point of view, on a grid of size h the smallest wavelength you can represent is 2*h and that requires to describe a sine by using at least three grid points.
FMDenaro is offline   Reply With Quote

Old   February 15, 2018, 17:11
Default
  #15
Senior Member
 
Join Date: Oct 2017
Location: United States
Posts: 233
Blog Entries: 1
Rep Power: 10
TurbJet is on a distinguished road
Quote:
Originally Posted by FMDenaro View Post
That's not exact. Your initial solution has still a spike in a point so that you cannot define the derivatives across it.
Also from a numerical point of view, on a grid of size h the smallest wavelength you can represent is 2*h and that requires to describe a sine by using at least three grid points.
I see your point.

I have been trying to manage all the terms, and here is what I get

LTE \approx -0.5 * \Delta t(a^{2}u_{xx} - abu_{xxx} + b^{2}u_{xxxx}) - \frac{a\Delta x^{2}}{6}u_{xxx} + \frac{b\Delta x^{2}}{12}u_{xxxx} + ...

but it seems the \Delta t only affect those terms in the brackets; can't see any terms in the form of \frac{\Delta x}{\Delta t}

I am confused again...
TurbJet is offline   Reply With Quote

Old   February 15, 2018, 17:14
Default
  #16
Senior Member
 
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 6,896
Rep Power: 73
FMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura about
I should check all the procedure, I remember some further terms...
FMDenaro is offline   Reply With Quote

Old   February 15, 2018, 17:22
Default
  #17
Senior Member
 
Join Date: Oct 2017
Location: United States
Posts: 233
Blog Entries: 1
Rep Power: 10
TurbJet is on a distinguished road
Quote:
Originally Posted by FMDenaro View Post
I should check all the procedure, I remember some further terms...
I attach my derivation below. Please check it if you have time.
TurbJet is offline   Reply With Quote

Old   February 15, 2018, 17:23
Default
  #18
Senior Member
 
Join Date: Oct 2017
Location: United States
Posts: 233
Blog Entries: 1
Rep Power: 10
TurbJet is on a distinguished road
Quote:
Originally Posted by FMDenaro View Post
I should check all the procedure, I remember some further terms...
Pic 1 derivation
Attached Images
File Type: jpg Screenshot from 2018-02-15 14:20:27.jpg (146.3 KB, 8 views)
TurbJet is offline   Reply With Quote

Old   February 15, 2018, 17:26
Default
  #19
Senior Member
 
Join Date: Oct 2017
Location: United States
Posts: 233
Blog Entries: 1
Rep Power: 10
TurbJet is on a distinguished road
Quote:
Originally Posted by FMDenaro View Post
I should check all the procedure, I remember some further terms...
pic 2 derivation
Attached Images
File Type: png Screenshot from 2018-02-15 14:23:50.png (95.7 KB, 10 views)
TurbJet is offline   Reply With Quote

Old   February 16, 2018, 13:25
Default
  #20
Senior Member
 
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 6,896
Rep Power: 73
FMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura about
I checked in my old notes and I extracted the LTE expression, it appears clearly that, provided that the dx allows to get the Reynolds cell number <=2, you have a global positive diffusion coefficient.
Now the issue is that if you compare the solution at some time T, when you decreases the time step you need more iteration to reach that time and the error effects summ.
Attached Files
File Type: pdf LTE di FTCS.pdf (148.7 KB, 5 views)
FMDenaro is offline   Reply With Quote

Reply

Tags
diffusion equation, time step size


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
pressure in incompressible solvers e.g. simpleFoam chrizzl OpenFOAM Running, Solving & CFD 13 March 28, 2017 06:49
simpleFoam error - "Floating point exception" mbcx4jc2 OpenFOAM Running, Solving & CFD 12 August 4, 2015 03:20
Help for the small implementation in turbulence model shipman OpenFOAM Programming & Development 25 March 19, 2014 11:08
Micro Scale Pore, icoFoam gooya_kabir OpenFOAM Running, Solving & CFD 2 November 2, 2013 14:58
How to write k and epsilon before the abnormal end xiuying OpenFOAM Running, Solving & CFD 8 August 27, 2013 16:33


All times are GMT -4. The time now is 19:07.