I have just tried the simple heat equation problem bellow:defined in the range x [0,1], t [0,1]with the initial condition: u(x,0)= 2x, 0<=x<=0.52(1-x), 0.5<x<=1.0 boundary condition:u(0,t)=u(1,t)=0I resort to Implicit Euler and Crank-Nicolson FDM respectivelythe results are shown in the following figure It seems that the Crank-Nicolson result is not very good.Comparing with Implicit Euler, C-N's result is not very smooth.I realize, Duffy's words:"Unfortunately, it has been known for some considerable time (Il’in, 1969) that centred differencing schemes in space combined with averaging in time (what essentially CN is in this context) leads to spurious oscillations in the approximate solution and in the divided differences for approximating its derivatives."Does this have anything to do with my results?Or just my results is wrong?thanks!

What a nice polished first post The IE is discussed here in two articles; they discuss why shows this behaviour and how IE with extrapolation achieves 2nd order accuracy and no oscillation.Extrapolation

Thank you!I will try your method.BTW, right now, I am reading your book and emailed you just a few day before for the code bug in the book

QuoteOriginally posted by: fogsnowThank you!I will try your method.BTW, right now, I am reading your book and emailed you just a few day before for the code bug in the book btw was that the post regarding the BAND matrix? Here is a thesis on the use of IE for Heston (scroll down to first post). http://www.datasimfinancial.com/forum/v ... .php?t=101

yesThanks for your sharing, I will check it.

You're welcome.On the related issue of matrix patterns, the boost uBLAS might be a good place to start.http://www.wilmott.com/blogs/cuchulainn/index.cfm

good suggestions! Next, I will learn boost and quantlib.

QuoteOriginally posted by: fogsnowI have just tried the simple heat equation problem bellow:defined in the range x [0,1], t [0,1]with the initial condition: u(x,0)= 2x, 0<=x<=0.52(1-x), 0.5<x<=1.0 boundary condition:u(0,t)=u(1,t)=0I resort to Implicit Euler and Crank-Nicolson FDM respectivelythe results are shown in the following figure It seems that the Crank-Nicolson result is not very good.Comparing with Implicit Euler, C-N's result is not very smooth.Here's a nice extra exercise: having computed U(j,n) for both IE and CN then calculate delta and gamma using divided difference, plot these values as above manner and watch what happens. CN is surprisingly beautiful Gibbs phenomenon. Most books do not show this pathological behaviour...

QuoteHere's a nice extra exercise: having computed U(j,n) for both IE and CN then calculate delta and gamma using divided difference, plot these values as above manner and watch what happens. CN is surprisingly beautiful Gibbs phenomenon. Most books do not show this pathological behaviour...

It is non-smooth initial condition (at x=0.5) that leads to the pathological behaviour.

QuoteOriginally posted by: fogsnowIt is non-smooth initial condition (at x=0.5) that leads to the pathological behaviour.Exactly. Just like a payoff function. CN looks like a fishess' spine, it's awful. We add it to the list of deprecated models.//a bit off, but I see that boost supports banded matrices. If you use my Matrix with structure then you could use boost band for it. I intend to do it this way and will put the code on my site.http://www.boost.org/doc/libs/1_38_0/li ... banded.htm

Quotea bit off, but I see that boost supports banded matrices. If you use my Matrix with structure then you could use boost band for it. I intend to do it this way and will put the code on my site.Ok, I will try boost.Your book is very helpful, the codes I used here are based on your book and I just make a little simplifications.

QuoteOriginally posted by: fogsnowIt is non-smooth initial condition (at x=0.5) that leads to the pathological behaviour.You are right. The roots of the problem you got with the CN scheme arethe non smooth (discoutinuous) initial data.This problem was solved in 1984 by RannacherRannacher, R (1984) Finite element solution of diffusion problems with irregulardata, Numerische Mathematik, 43, 309-327.Simply replace the CN approximation for the first timestep by two half-timesteps usingbackward Euler time integration, eg(I + 1/2 K_j) U_{j}^{n+1} = (I - 1/2K_j) U_{j}^{n}by(I + 1/2 K_j) U_{j}^{n+1} = U_{j}^{n}where K is the discretisation of the diffusive, convective... terms.This means taking two fully implicit methods at the start of CN. This leads to an oscillatons-free scheme.You can also have a look to this paperR. Carter and M.B. Giles. 'Sharp error estimates for discretisations of the 1D convection/diffusion equation with Dirac initial data', IMA Journal of Numerical Analysis, 27(2):406-425, 2007Farid.PS: I already posted this "answer" on Monday Dec 01, 08 10:52 AM in the thread "crank nicholson code", Numerical Methods Forum.

Unfortunately, Rannacher method is not robust under all parameter values, for sure.Why not just accept extrapolation and be done with it??Now Rannacher == IE for the first few steps, then CN after that. This is a kind of workaround as far as I am concerned. QuoteThe roots of the problem you got with the CN scheme are the non smooth (discoutinuous) initial data.Yes, this is being hammered to death at this stage. //BTW"crank nicholson code", should be "nicolson"

QuoteThis means taking two fully implicit methods at the start of CN. This leads to an oscillatons-free scheme.This is scary1. how many steps at the beginning to avoid OSC?2. the consequences of going from first order to second-order accuracy is ad-hoc in my opinion.3. extra complicated code4. it's a hassle5. von Neumann stable < L0 stable (can handle discontinuous IC)You can read Sheppards' (Heston) thesis on Datasim site and why he did NOT use CN. Don't just take my word for it