C++ 3D PDE FDM

[xR-]

Dear all,

I am currently trying to solve high-dimension PDEs for a school project.
The Boost and Computational Finance page is very helpful for me, especially the 4th part !

I tried to gather the different sections of code in Visual Studio, to add the missing #include, as well as the snippets from Duffy's 2006 CD-Rom, but I still have a few issues ! It will build, but crash at run time, probably because I did not initialized correctly a few BoostTensors.

My whole project is about pricing with a rates/equity hybrid diffusion, unfortunately I can not afford spending months to study in details the different PDE solving methods but the one introduced in the article really sounds great so I would be more than happy to use this one !

I think I am really close to get it working but I can not finish :-/ Do someone here (Daniel Duffy) ? would have by any chance a running code of this solver + heat example ?

Many thanks in advance ,

Very best.

 

[Daniel Duffy]

What pde do you wish to solve, exactly?

 

[xR-]

Kolmogorov Forward on a 2D density

(an alternative 2D method could be welcome, but it would be nice to keep the 3D feature for future possible evolutions)

EDIT> Yes, maybe the high-dimension mention in my 1st post was a bit over-emphasized ;-)

 

[Daniel Duffy]

Ah, so 2 space variables and one time.

What's next?

e.g. Writing up the pde is a good start.

 

[xR-]

Hello Daniel,

The exact PDE is the following :

(\frac{\partial \psi_t}{\partial t}+ (r_t-f(0,t))\psi_t + \frac{\partial (r_t S_t \psi_t)}{\partial S} + \frac{\partial (a_t \psi_t)}{\partial r} -\frac{1}{2} \frac{\partial^2 (\sigma_t^2 S_t^2 \psi_t)}{\partial S^2} -\frac{1}{2} \frac{\partial^2 (b_t^2 \psi_t)}{\partial r^2} - \rho \frac{\partial^2 (\sigma_t b_t S_t \psi_t)}{\partial S \partial r} = 0)

where :

a, b denotes short rate diffusion coefficients, rho is the constant correlation between two 1d-wiener processes and psi is the 2d-density that we try to propagate in time.

sigma is the equity volatility, we need to solve it recursively, alternatively with the 2d-density (by computing numerically an expectation, with respect to this 2d-density)

We might need the 3d if we switch to a 2 factors rates diffusion (for now we use a 1d Hull-White))

Regarding the code published on QuantNet, I have some difficulties because of the solver class : the initialization of some members + writing correctly the loop over time :-/ I am aware it is not ultra difficult, but it would greatly help if a self-contained source-code was available, as a starting point.

Best,

 

[Daniel Duffy]

xR-,
This is the most ambitious post that I have ever seen on Quantnet I am having some difficulty with the scope of your project.

Some questions:

1. What is your goal with this endeavour?
2. Do you have any exposure to PDE/FDM?
3. The code on QN avoids many essential complexities that your PDE has (too many to enumerate here)
4. How long is this project and what is your contribution and what do you wish to use from others?

Here is a thesis on Heston to give an indication of what is required. (go to first post where you can download thesis).

http://www.datasimfinancial.com/forum/viewtopic.php?t=101

See also FD2 prioject here which is not a million miles away from your current PDE because it used ADE.

http://qfcl.wilmott.com/

Let me know what you think.

Thanks

Daniel

 

[DStahl]

I am trying to do something similar (though only in one dimension...my goal is to do more efficient Monte Carlo simulations by using the numerical distribution instead of simulating the path of the process). While I have a numerical solution that "looks" good, the area under the curve at any point in time does not sum to one. In fact it doesn't even seem to converge to one as the time and space steps get very small. Perhaps even worse, the area under the curve varies for different points in time. I wonder how people get around this?

 

[Daniel Duffy]

I am trying to do something similar (though only in one dimension...my goal is to do more efficient Monte Carlo simulations by using the numerical distribution instead of simulating the path of the process). While I have a numerical solution that "looks" good, the area under the curve at any point in time does not sum to one. In fact it doesn't even seem to converge to one as the time and space steps get very small. Perhaps even worse, the area under the curve varies for different points in time. I wonder how people get around this?

You are computing the transition density and then use in numerical quadrature?
Which FDM did you use?

 

[DStahl]

I used implicit on the GBM stock price, since a closed form solution is available and I can check my answers. I used a normal density with variance sigma^2 delta t, and mean x+alpha delta t as my initial condition since I was getting weird values when using a single point as my IC. I then did numerical quadrature on various points on the time axis. This is what the surface plot looks like, which seems about right.

 

[DStahl]

I finally got it to work. For some reason the forward equation doesn't seem to work if I take the derivatives before discretizing the PDE, but it does if I discretize the original PDE. Or maybe I just made a mistake taking the derivatives (though that seems unlikely, I triple checked my work and the partials are not hard to solve).

 

[Daniel Duffy]

I finally got it to work. For some reason the forward equation doesn't seem to work if I take the derivatives before discretizing the PDE, but it does if I discretize the original PDE. Or maybe I just made a mistake taking the derivatives (though that seems unlikely, I triple checked my work and the partials are not hard to solve).

Most (all?) QF articles take the derivative of (ap)_ss + (bp)_s + .. because most fdm is Crank Nicolson.

The reason is usually singularities in diffusion etc.

We can also keep the original PDE and apply ADE method, for example. It means no LU decomposition needed at each level.

BTW do you have a modified tridiagonal matrix?

 

[DStahl]

What do you mean modified? It contains the coefficients of the difference equation, which are space dependent in this case. Since the space boundaries are equal to zero I haven't done anything for those BCs.

 

[Daniel Duffy]

What do you mean modified? It contains the coefficients of the difference equation, which are space dependent in this case. Since the space boundaries are equal to zero I haven't done anything for those BCs.

let v = ap, Then the diffusion becomes

(v(j+1) - 2v(j) + v(j-1)/h^2 --> (a(j+1)U(j+1) ... etc.)

I think we are saying the same thing; 'modified' was not a good word.