I'm a beginner with FE. My application is the pricing of financial derivatives where the space is five dimensional. So, adding time, the problem has six dimensions.

I tried to look around (Fenics, escript, deal.II,...), but my understanding is that those software are limited to 3 + 1 (3d space + 1d time). Is this correct?

My targeted language are Python or C++.

Description of my problem
I would like to price an investment product where, each month, the investor has the freedom to re-invests or not. I would like to do so with stochastic volatility, stochastic interest rate, and stochastic mortality.
The stochastic PDEs look like this \begin{align} dS_t &= \mu^S_t d_t + \sqrt{\sigma_t} dB^S_t &\text{(stock)}\\ d\sigma_t &= \mu^\sigma_t dt + \nu^\sigma_t dB^\sigma_t & \text{(volatility)} \\ dr_t &= \mu^r_t dt + \nu^r_t dB^r_t & \text{(interest rate)} \\ dq_t &= \mu^q_t dt + \nu^q_t dB^q_t & \text{(mortality)} \end{align} Where $\mu^S_t$ is a time-dependent constant associated to the stock price $S$, and $B^S_t$ is an independent Levy process which creates noise in the stock price $S$. Similarly for the other quantities: $\nu^\sigma_t$ is a time-dependent quantity associated to the volatility $\sigma$.
Let $C_\tau$ denotes the admissible investments at time $\tau$. The stochastic control problem looks like $$ V_\tau = max \left\{ c \in C_\tau : P(\text{death})E(r_\tau f(S_{\tau+1})) + P(alive)E(r_\tau V_{\tau+1})\right\}. $$ The above PDEs are continuous, but the value of the product $V_\tau$ is solved only at predefined $\tau$-times, say each month.

I guess Monte-Carlo can always brute force my problem, but it is very slow.

Deterministic form of the stochastic PDEs
For this part, assume that the value of the option $$ V : (t, S_t, \sigma_t, r_t, q_t, c_t) \mapsto (t, V_t), $$ is define on the natural time $t$, not the $\tau$-times, with $c_t$ the investment at time $t$.
Define the differential operator \begin{align} L_t &= \partial_{r,S} + \partial_{r,\sigma} + \partial_{\sigma,S} \\ L^S_t &= \sigma_t \partial_S + r_t \partial_{S,S} \\ L^r_t &= \partial_r + \partial_{r,r} \\ L^\sigma_t &= \partial_\sigma + \partial_{\sigma,\sigma} \\ L^q_t &= \partial_q + \partial_{q,q} \end{align} where time-dependent constant $\{\mu^S_t,\ldots\}$ are ignored. The deterministic PDE is then $$ \partial_t V_t +\left(L_t+ L^S_t + L^\sigma_t + L^r_t+L^q_t\right)V_t = 0, $$ which can adapted to the optimal control problem on the $\tau$-times.

  • 3
    $\begingroup$ Are you sure you need to use finite elements for this problem? It would help if you can describe the problem some more (in particular the PDE that you want to solve). $\endgroup$
    – Victor Liu
    Commented Sep 10, 2013 at 23:51
  • $\begingroup$ @Liu I added more details. I though about FE because MC is very slow. $\endgroup$
    – user729
    Commented Sep 11, 2013 at 0:31
  • 1
    $\begingroup$ Can you clarify the notation? Does $v^p$ imply a derivative in $p$? $\endgroup$
    – Jesse Chan
    Commented Sep 11, 2013 at 0:47
  • 3
    $\begingroup$ I think you'll get better answers if you also post the deterministic PDEs you're going to solve. Could you clarify what the independent variables are? Right now, it looks like the only independent variable is time. Are you solving these stochastic differential equations using polynomial chaos expansions, and is that why you will have a system of deterministic differential equations? $\endgroup$ Commented Sep 11, 2013 at 6:05
  • 1
    $\begingroup$ On the one hand you might deal with the complications of using FEs in moderate dimensions and the curse of dimensionality, or you can work on speedup methods for MC or better QMC. The latter world is not necessarily worse, actually it is the approach of choice in the quant world for many reasons, so be careful in dismissing it so easily. $\endgroup$
    – Quartz
    Commented Sep 13, 2013 at 9:23

3 Answers 3


Assuming you want to solve the Black-Scholes equations or a variant on a portfolio of 5 assets, then you indeed have 5 spatial plus one time dimension. I don't know any FEM package that can do that off the top of my head (deal.II can't readily do it, but see below) but I think I remember that some people from Chris Schwab's group at ETH Zurich solved such problems using sparse meshes. You may get lucky looking around his publications.

There are other equations that have extra dimensions. One example is the radiative transfer equation that has 3 space + 1 time + 2 angular + 1 energy dimension. The way this is typically solved is to discretize 3-dimensional space as usual, then discretize the angular and energy dimensions on separate 2 and 1-dimensional meshes and at each nodal point of the spatial mesh simply have a lot of variables (one each for each node of the angular mesh times the number of nodes in the energy mesh). We use this scheme in deal.II implementations successfully. This makes sense for the radiative transfer equation, and it may be emulated for your equation even if it isn't natural there.


DUNE, the Distributed and Unified Numerics Environment http://www.dune-project.org, features some structured grids of arbitrary dimension (SGrid, and Yaspgrid), see DUNE's features. Currently, there is a branch that turns yaspgrid, one of the above grids, into a tensor product grid, if that is of interest. Since release 2.0 (the current one is 2.2.1, and 2.3 is soon to come) we do have reference elements for various finite element methods that support arbitrary dimensions. Therefore it should be possible set up finite element discretizations of arbitrary dimension with e.g. the disrectization module dune-pdelab. Albeit this might not be tested to often.

Having said that, there is still the curse of dimensionality as Wolfgang pointed out.

For further information I refer you to the DUNE mailinglists.


Ok, so it looks like what you have are a coupled set of ODEs, since as far as I can tell, there are only derivatives with respect to time, and no derivatives with respect to anything else. There are quite a few packages out there for solving systems of ODEs of arbitrary dimension (Matlab has stuff like ode45). For Python, look at this question for some suggestions. Finally, there is old Fortran code on netlib that can be interfaced with C++ pretty easily (ease of use is another matter). There are probably better alternatives out there since it's been a while since I've looked (others should chime in).

  • 1
    $\begingroup$ By adding the deterministic PDEs, I see my question was not clear. Sorry and thank you for trying to help. $\endgroup$
    – user729
    Commented Sep 11, 2013 at 14:38

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.