Solving Dirichlet problems numerically using the Feynman-Kac representation

Fabian M. Buchmann, Wesley P. Petersen

Seminar für Angewandte Mathematik, ETH Zürich, ETH Zentrum, CH-8092 Zürich, Switzerland

fab@math.ethz.ch

Contributed talk

We study numerical solutions of the Dirichlet problem in high dimensions using the Feynman-Kac representation. What is involved are Monte-Carlo simulations of stochastic differential equations and algorithms to accurately determine exit times and process values at the boundary. In a basic form, the Feynman-Kac representation for the Dirichlet problem is follows. Let $L$ be a linear second order differential operator and consider the PDE

$$Lu = \frac{1}{2}\bigtriangleup u+\sum_{i=1}^n b_i(x)\partial_iu +c(x)u=-g(x).$$ (1)

We assume that the potential $c(x)\le 0$ and that $b,c,g$ are smooth and satisfy some Lipshitz growth conditions. The Dirichlet Problem is to find $u(x)$ satisfying above PDE (1), while at the boundary of a bounded domain $D$, $u(x)=\psi(x)$, when $x\in\partial D$. The probabilistic representation for the solution $u$ is the Feynman-Kac formula,

$$u(x)=\mathbf{E}\psi(X(\tau))e^{\int_0^\tau c(X(s))ds} + \mathbf{E}\int_0^\tau g(X(t))e^{\int_0^tc(X(s))}ds dt,$$ (2)

where $\tau$ is the exit time from domain $D$ of process $X(t)$ whose initial value is $X(0)=x$ and which satisfies a stochastic differential equation $dX(t)=b(X)dt+dW$.

Two difficulties present themselves when evaluating (2): (i) finding the process value $X(\tau)$ accurately at the boundary, and (ii) making sure that this value does not follow an excursion. If the numerical representation of the increments of $W$ is Gaussian, controlling both problems is very hard because these increments are unbounded. We show how to solve this problems using a 3-point approximation to $\Delta W$ and explain the resulting walk on cubes. We study problems in a hypersherical domain on the basic test problem

$$\frac{1}{2}\bigtriangleup u + g(x) = 0, \quad u(x)\big\vert _{\partial D}=\psi(x)$$

in high dimensions (with $n$ up to 64). Our particular goal is to see, how the solutions behave as the dimension $n$ of the problem increases. We compare our approach with other solutions, namely a walk on spheres.

We find that the canonical $O(N^{-1/2})$ behavior of statistical errors as a function of the sample size $N$ holds regardless of the dimension $n$ of the space. In fact, the coefficient of $N^{-1/2}$ seems to actually decrease with $n$.