Namespaces
Variants
Actions

Difference between revisions of "Monte-Carlo methods for partial differential equations"

From Encyclopedia of Mathematics
Jump to: navigation, search
(Importing text file)
 
m (tex encoded by computer)
 
Line 1: Line 1:
 +
<!--
 +
m1101901.png
 +
$#A+1 = 66 n = 0
 +
$#C+1 = 66 : ~/encyclopedia/old_files/data/M110/M.1100190 Monte\AAhCarlo methods for partial differential equations
 +
Automatically converted into TeX, above some diagnostics.
 +
Please remove this comment and the {{TEX|auto}} line below,
 +
if TeX found to be correct.
 +
-->
 +
 +
{{TEX|auto}}
 +
{{TEX|done}}
 +
 
Deterministic partial differential equations can be solved numerically by probabilistic algorithms such as Monte-Carlo methods, stochastic particle methods, ergodic algorithms, etc. (cf. also [[Differential equation, partial|Differential equation, partial]]).
 
Deterministic partial differential equations can be solved numerically by probabilistic algorithms such as Monte-Carlo methods, stochastic particle methods, ergodic algorithms, etc. (cf. also [[Differential equation, partial|Differential equation, partial]]).
  
The simplest example is the [[Heat equation|heat equation]] in <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m1101901.png" />:
+
The simplest example is the [[Heat equation|heat equation]] in $  C ^ {1,2 } ( ( 0,T ] \times \mathbf R  ^ {d} ) $:
  
<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m1101902.png" /></td> </tr></table>
+
$$
 +
\left \{
 +
\begin{array}{l}
 +
{ {
 +
\frac{\partial  u }{\partial  t }
 +
} = {
 +
\frac{1}{2}
 +
} \Delta u, \  } \\
 +
{u ( t,x ) \rightarrow u _ {0} ( x ) \  \textrm{ as  }  t \rightarrow 0, } \\
 +
{\  \textrm{ for  all  }  x  \textrm{ at  which  } \ , u _ {0}  \textrm{ is  continuous  } . }
 +
\end{array}
 +
\right .
 +
$$
  
For a large class of functions <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m1101903.png" />, the unique solution is <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m1101904.png" />, where <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m1101905.png" /> is the heat kernel. Thus, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m1101906.png" />, where <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m1101907.png" /> is any [[Brownian motion|Brownian motion]] such that <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m1101908.png" />, almost surely. The strong [[Law of large numbers|law of large numbers]] applies: if <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m1101909.png" /> are independent copies of <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019010.png" />, one has, almost surely,
+
For a large class of functions $  u _ {0} $,  
 +
the unique solution is $  u ( t,x ) = g _ {t} * u _ {0} ( x ) $,  
 +
where $  g _ {t} $
 +
is the heat kernel. Thus, $  u ( t,x ) = {\mathsf E} u _ {0} ( W _ {t} ( x ) ) $,  
 +
where $  ( W _ {t} ( x ) ) $
 +
is any [[Brownian motion|Brownian motion]] such that $  W _ {0} = x $,  
 +
almost surely. The strong [[Law of large numbers|law of large numbers]] applies: if $  \{ ( W _ {t}  ^ {i} ( x ) ) \} $
 +
are independent copies of $  ( W _ {t} ( x ) ) $,  
 +
one has, almost surely,
  
<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019011.png" /></td> </tr></table>
+
$$
 +
u ( t,x ) = {\lim\limits } _ {N \rightarrow \infty } {
 +
\frac{1}{N}
 +
} \sum _ {i = 1 } ^ { N }  u _ {0} ( W _ {t}  ^ {i} ( x ) ) .
 +
$$
  
A Monte-Carlo procedure to compute <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019012.png" /> at a single point <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019013.png" /> consists in the simulation of <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019014.png" /> independent copies of <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019015.png" /> and the computation of
+
A Monte-Carlo procedure to compute $  u ( t, \cdot ) $
 +
at a single point $  x $
 +
consists in the simulation of $  N $
 +
independent copies of $  W _ {t} ( x ) $
 +
and the computation of
  
<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019016.png" /></td> </tr></table>
+
$$
 +
{
 +
\frac{1}{N}
 +
} \sum _ {i = 1 } ^ { N }  u _ {0} ( W _ {t}  ^ {i} ( x ) ) .
 +
$$
  
This procedure requires the simulation of <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019017.png" /> independent real-valued Gaussian random variables of mean <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019018.png" /> and variance <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019019.png" /> (cf. [[Random variable|Random variable]]; [[Gaussian process|Gaussian process]]). The convergence rate of the procedure is described by fluctuation results, large-deviation estimates, etc. (cf. also [[Probability of large deviations|Probability of large deviations]]). The most elementary estimate concerns the standard deviation of the random error:
+
This procedure requires the simulation of $  N \times d $
 +
independent real-valued Gaussian random variables of mean 0 $
 +
and variance $  t $(
 +
cf. [[Random variable|Random variable]]; [[Gaussian process|Gaussian process]]). The convergence rate of the procedure is described by fluctuation results, large-deviation estimates, etc. (cf. also [[Probability of large deviations|Probability of large deviations]]). The most elementary estimate concerns the standard deviation of the random error:
  
<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019020.png" /></td> </tr></table>
+
$$
 +
{
 +
\frac{1}{\sqrt N }
 +
} \sqrt { {\mathsf E} \left | {u _ {0} ( W _ {t} ( x ) ) } \right |  ^ {2} - \left | {u ( t,x ) } \right |  ^ {2} } .
 +
$$
  
More generally, let <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019021.png" /> be an integro-differential operator (cf. also [[Integro-differential equation|Integro-differential equation]]), and consider the evolution problem
+
More generally, let $  L $
 +
be an integro-differential operator (cf. also [[Integro-differential equation|Integro-differential equation]]), and consider the evolution problem
  
<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019022.png" /></td> </tr></table>
+
$$
 +
{
 +
\frac{d}{dt }
 +
} u ( t,x ) = Lu ( t,x ) ,  u ( 0,x ) = u _ {0} ( x ) .
 +
$$
  
From now on, it is assumed that there is a <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019023.png" />-valued [[Markov process|Markov process]] <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019024.png" /> such that <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019025.png" /> (for a discussion on the assumptions that <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019026.png" /> must satisfy, see [[#References|[a3]]], Chap. XIII, Sect. 3). Then, if it is possible to simulate independent copies of the process <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019027.png" />, one can approximate <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019028.png" /> by
+
From now on, it is assumed that there is a $  \mathbf R  ^ {d} $-
 +
valued [[Markov process|Markov process]] $  ( X _ {t} ) $
 +
such that $  u ( t,x ) = {\mathsf E} u ( 0,X _ {t} ( x ) ) $(
 +
for a discussion on the assumptions that $  L $
 +
must satisfy, see [[#References|[a3]]], Chap. XIII, Sect. 3). Then, if it is possible to simulate independent copies of the process $  ( X _ {t} ( x ) ) $,  
 +
one can approximate $  u ( t,x ) $
 +
by
  
<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019029.png" /></td> </tr></table>
+
$$
 +
{
 +
\frac{1}{N}
 +
} \sum _ {i = 1 } ^ { N }  u _ {0} ( X  ^ {i} _ {t} ( x ) ) .
 +
$$
  
For example, in neutronics, the process <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019030.png" /> is the pair (position,velocity), and the velocity is a pure [[Jump process|jump process]] (in this case, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019031.png" /> is a transport operator) whose simulation is simple: one simulates the times at which the velocity changes and the new velocity after each change time independently (the change times have an exponential law).
+
For example, in neutronics, the process $  ( X _ {t} ) $
 +
is the pair (position,velocity), and the velocity is a pure [[Jump process|jump process]] (in this case, $  L $
 +
is a transport operator) whose simulation is simple: one simulates the times at which the velocity changes and the new velocity after each change time independently (the change times have an exponential law).
  
Most often, the law of <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019032.png" /> cannot be simulated exactly, but it may be possible to approximate it suitably. For example, when <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019033.png" /> is the solution of a [[Stochastic differential equation|stochastic differential equation]] driven by a Brownian motion (in this case, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019034.png" /> is a second-order differential operator), one can discretize the stochastic differential equation in time. The simplest discretization scheme (the Euler scheme) defines an approximate process <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019035.png" /> depending on the discretization step <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019036.png" />, whose simulation requires the simulation of independent Gaussian random variables only. The approximation formula is
+
Most often, the law of $  X _ {t} ( x ) $
 +
cannot be simulated exactly, but it may be possible to approximate it suitably. For example, when $  ( X _ {t} ( x ) ) $
 +
is the solution of a [[Stochastic differential equation|stochastic differential equation]] driven by a Brownian motion (in this case, $  L $
 +
is a second-order differential operator), one can discretize the stochastic differential equation in time. The simplest discretization scheme (the Euler scheme) defines an approximate process $  ( X  ^ {h} _ {t} ( x ) ) $
 +
depending on the discretization step $  h $,  
 +
whose simulation requires the simulation of independent Gaussian random variables only. The approximation formula is
  
<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019037.png" /></td> </tr></table>
+
$$
 +
u ( t,x ) \simeq {
 +
\frac{1}{N}
 +
} \sum _ {i = 1 } ^ { N }  u _ {0} ( X _ {t} ^ {h,i } ( x ) ) .
 +
$$
  
 
The global error of this procedure can be decomposed as the sum of
 
The global error of this procedure can be decomposed as the sum of
  
<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019038.png" /></td> </tr></table>
+
$$
 +
u ( t,x ) - {\mathsf E} u _ {0} ( X  ^ {h} _ {t} ( x ) )
 +
$$
  
 
(the discretization error) and of
 
(the discretization error) and of
  
<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019039.png" /></td> </tr></table>
+
$$
 +
{\mathsf E} u _ {0} ( X  ^ {h} _ {t} ( x ) ) - {
 +
\frac{1}{N}
 +
} \sum _ {i = 1 } ^ { N }  u _ {0} ( X _ {t} ^ {h,i } ( x ) )
 +
$$
  
(the statistical error). Under appropriate assumptions, the discretization error can be expanded in terms of <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019040.png" />.
+
(the statistical error). Under appropriate assumptions, the discretization error can be expanded in terms of $  h $.
  
 
The numerical procedure can be adapted to the solution of problems in bounded domains. Dirichlet boundary conditions lead to the simulation of processes stopped at the boundary. The simulation of reflected diffusion processes (cf. [[Diffusion process|Diffusion process]]) solves parabolic or elliptic partial differential equations with suitable [[Neumann boundary conditions|Neumann boundary conditions]].
 
The numerical procedure can be adapted to the solution of problems in bounded domains. Dirichlet boundary conditions lead to the simulation of processes stopped at the boundary. The simulation of reflected diffusion processes (cf. [[Diffusion process|Diffusion process]]) solves parabolic or elliptic partial differential equations with suitable [[Neumann boundary conditions|Neumann boundary conditions]].
  
Let <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019041.png" /> be the formal adjoint of <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019042.png" />. Note that the probability distribution <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019043.png" /> of <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019044.png" /> is a solution to
+
Let $  L  ^ {*} $
 +
be the formal adjoint of $  L $.  
 +
Note that the probability distribution $  \mu _ {t}  ^ {x} $
 +
of $  X _ {t} ( x ) $
 +
is a solution to
  
<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019045.png" /></td> </tr></table>
+
$$
 +
{
 +
\frac{d}{dt }
 +
} \mu _ {t}  ^ {x} = L  ^ {*} \mu _ {t}  ^ {x}
 +
$$
  
in the sense of distributions (cf. [[Generalized function|Generalized function]]). The solution <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019046.png" /> can be approximated by the [[Empirical distribution|empirical distribution]]
+
in the sense of distributions (cf. [[Generalized function|Generalized function]]). The solution $  \mu _ {t}  ^ {x} $
 +
can be approximated by the [[Empirical distribution|empirical distribution]]
  
<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019047.png" /></td> </tr></table>
+
$$
 +
{
 +
\frac{1}{N}
 +
} \sum _ {i = 1 } ^ { N }  \delta _ {X _ {t}  ^ {i} ( x ) } .
 +
$$
  
 
This procedure is used in neutronics.
 
This procedure is used in neutronics.
  
To certain non-linear partial differential equations one can associate processes <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019048.png" /> which are not Markovian, but are such that the pairs <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019049.png" /> are Markovian, where <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019050.png" /> is the probability law of the random variable <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019051.png" />. Besides, in some cases, the law <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019052.png" /> of the process <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019053.png" /> can be described by means of the theory of propagation of chaos (see [[#References|[a4]]] and the references therein). That means that the law <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019054.png" /> is the limit in probability, when <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019055.png" /> goes to infinity, of the empirical distributions of the paths of <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019056.png" /> particles in a mean field interaction. The interaction between the particles is related to the non-linear differential operator of the partial differential equation under consideration. If one can simulate the stochastic particle system, then an algorithm can be built from the formula
+
To certain non-linear partial differential equations one can associate processes $  ( X _ {t} ) $
 +
which are not Markovian, but are such that the pairs $  ( X _ {t} , \mu _ {t} ) $
 +
are Markovian, where $  \mu _ {t} $
 +
is the probability law of the random variable $  X _ {t} $.  
 +
Besides, in some cases, the law $  P  ^ {X} $
 +
of the process $  ( X _ {t} ) $
 +
can be described by means of the theory of propagation of chaos (see [[#References|[a4]]] and the references therein). That means that the law $  P  ^ {X} $
 +
is the limit in probability, when $  N $
 +
goes to infinity, of the empirical distributions of the paths of $  N $
 +
particles in a mean field interaction. The interaction between the particles is related to the non-linear differential operator of the partial differential equation under consideration. If one can simulate the stochastic particle system, then an algorithm can be built from the formula
  
<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019057.png" /></td> </tr></table>
+
$$
 +
\mu _ {t} \simeq {
 +
\frac{1}{N}
 +
} \sum _ {i = 1 } ^ { N }  \delta _ {X _ {t}  ^ {i} } ,
 +
$$
  
where <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019058.png" /> are the locations at time <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019059.png" /> of the interacting particles. For example, consider an interaction kernel <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019060.png" />, a family of <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019061.png" /> independent Brownian motions, and the system
+
where $  \{ X  ^ {i} _ {t} \} $
 +
are the locations at time $  t $
 +
of the interacting particles. For example, consider an interaction kernel $  K $,  
 +
a family of $  N $
 +
independent Brownian motions, and the system
  
<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019062.png" /></td> </tr></table>
+
$$
 +
X _ {t}  ^ {i} = X _ {0}  ^ {i} + \int\limits _ { 0 } ^ { t }  { {
 +
\frac{1}{N}
 +
} \sum _ {j = 1 } ^ { N }  K ( X _ {s}  ^ {i} - X _ {s}  ^ {j} ) }  {ds } + \nu W _ {t}  ^ {i} .
 +
$$
  
The time discretization of such systems permits one to approximate equations such as Vlasov–Poisson–Fokker–Planck equations (cf. [[Vlasov–Poisson–Fokker–Planck system|Vlasov–Poisson–Fokker–Planck system]]), the viscous Burgers equation (cf. [[Diffusion equation|Diffusion equation]]; [[Turbulence, mathematical problems in|Turbulence, mathematical problems in]]), the incompressible Navier–Stokes equation for the vorticity in <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019063.png" /> (cf. [[Navier–Stokes equations|Navier–Stokes equations]]), etc. In the latter case, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019064.png" /> is the Biot–Savart kernel, and the algorithm is the Chorin random vortex method (see [[#References|[a2]]] for numerical considerations). Similarly, the simulation of interacting branching processes (cf. [[Branching process|Branching process]]) permits one to approximate the solutions of convection-reaction-diffusion equations. One difficulty of the mathematical study of these methods is to show that the statistical errors (on the distribution function of <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019065.png" />, e.g.) have standard deviations of order <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/m/m110/m110190/m11019066.png" /> although the particles are dependent. For theoretical aspects related to probabilistic algorithms for partial differential equations, including the [[Boltzmann equation|Boltzmann equation]], see, for instance, [[#References|[a1]]] and the references therein.
+
The time discretization of such systems permits one to approximate equations such as Vlasov–Poisson–Fokker–Planck equations (cf. [[Vlasov–Poisson–Fokker–Planck system|Vlasov–Poisson–Fokker–Planck system]]), the viscous Burgers equation (cf. [[Diffusion equation|Diffusion equation]]; [[Turbulence, mathematical problems in|Turbulence, mathematical problems in]]), the incompressible Navier–Stokes equation for the vorticity in $  \mathbf R  ^ {2} $(
 +
cf. [[Navier–Stokes equations|Navier–Stokes equations]]), etc. In the latter case, $  K $
 +
is the Biot–Savart kernel, and the algorithm is the Chorin random vortex method (see [[#References|[a2]]] for numerical considerations). Similarly, the simulation of interacting branching processes (cf. [[Branching process|Branching process]]) permits one to approximate the solutions of convection-reaction-diffusion equations. One difficulty of the mathematical study of these methods is to show that the statistical errors (on the distribution function of $  \mu _ {t} $,  
 +
e.g.) have standard deviations of order $  {1 / {\sqrt N } } $
 +
although the particles are dependent. For theoretical aspects related to probabilistic algorithms for partial differential equations, including the [[Boltzmann equation|Boltzmann equation]], see, for instance, [[#References|[a1]]] and the references therein.
  
 
====References====
 
====References====
 
<table><TR><TD valign="top">[a1]</TD> <TD valign="top">  C. Graham,  T. Kurtz,  S. Méléard,  P. Protter,  M. Pulvirenti,  D. Talay,  "Probabilistic models for nonlinear partial differential equations" , ''Lecture Notes in Mathematics'' , '''1627''' , Springer  (1996)</TD></TR><TR><TD valign="top">[a2]</TD> <TD valign="top">  "Vortex methods and vortex motions"  K.E. Gustavson (ed.)  J.A. Sethian (ed.) , SIAM  (1991)</TD></TR><TR><TD valign="top">[a3]</TD> <TD valign="top">  J. Jacod,  "Calcul stochastique et problèmes de martingales" , ''Lecture Notes in Mathematics'' , '''714''' , Springer  (1979)</TD></TR><TR><TD valign="top">[a4]</TD> <TD valign="top">  A.S. Sznitman,  "Topics in propagation of chaos"  P.L. Hennequin (ed.) , ''Ecole d'Eté de Probabilités de Saint-Flour XI (1989)'' , ''Lecture Notes in Mathematics'' , '''1464''' , Springer  (1991)  pp. 165–251</TD></TR></table>
 
<table><TR><TD valign="top">[a1]</TD> <TD valign="top">  C. Graham,  T. Kurtz,  S. Méléard,  P. Protter,  M. Pulvirenti,  D. Talay,  "Probabilistic models for nonlinear partial differential equations" , ''Lecture Notes in Mathematics'' , '''1627''' , Springer  (1996)</TD></TR><TR><TD valign="top">[a2]</TD> <TD valign="top">  "Vortex methods and vortex motions"  K.E. Gustavson (ed.)  J.A. Sethian (ed.) , SIAM  (1991)</TD></TR><TR><TD valign="top">[a3]</TD> <TD valign="top">  J. Jacod,  "Calcul stochastique et problèmes de martingales" , ''Lecture Notes in Mathematics'' , '''714''' , Springer  (1979)</TD></TR><TR><TD valign="top">[a4]</TD> <TD valign="top">  A.S. Sznitman,  "Topics in propagation of chaos"  P.L. Hennequin (ed.) , ''Ecole d'Eté de Probabilités de Saint-Flour XI (1989)'' , ''Lecture Notes in Mathematics'' , '''1464''' , Springer  (1991)  pp. 165–251</TD></TR></table>

Latest revision as of 08:01, 6 June 2020


Deterministic partial differential equations can be solved numerically by probabilistic algorithms such as Monte-Carlo methods, stochastic particle methods, ergodic algorithms, etc. (cf. also Differential equation, partial).

The simplest example is the heat equation in $ C ^ {1,2 } ( ( 0,T ] \times \mathbf R ^ {d} ) $:

$$ \left \{ \begin{array}{l} { { \frac{\partial u }{\partial t } } = { \frac{1}{2} } \Delta u, \ } \\ {u ( t,x ) \rightarrow u _ {0} ( x ) \ \textrm{ as } t \rightarrow 0, } \\ {\ \textrm{ for all } x \textrm{ at which } \ , u _ {0} \textrm{ is continuous } . } \end{array} \right . $$

For a large class of functions $ u _ {0} $, the unique solution is $ u ( t,x ) = g _ {t} * u _ {0} ( x ) $, where $ g _ {t} $ is the heat kernel. Thus, $ u ( t,x ) = {\mathsf E} u _ {0} ( W _ {t} ( x ) ) $, where $ ( W _ {t} ( x ) ) $ is any Brownian motion such that $ W _ {0} = x $, almost surely. The strong law of large numbers applies: if $ \{ ( W _ {t} ^ {i} ( x ) ) \} $ are independent copies of $ ( W _ {t} ( x ) ) $, one has, almost surely,

$$ u ( t,x ) = {\lim\limits } _ {N \rightarrow \infty } { \frac{1}{N} } \sum _ {i = 1 } ^ { N } u _ {0} ( W _ {t} ^ {i} ( x ) ) . $$

A Monte-Carlo procedure to compute $ u ( t, \cdot ) $ at a single point $ x $ consists in the simulation of $ N $ independent copies of $ W _ {t} ( x ) $ and the computation of

$$ { \frac{1}{N} } \sum _ {i = 1 } ^ { N } u _ {0} ( W _ {t} ^ {i} ( x ) ) . $$

This procedure requires the simulation of $ N \times d $ independent real-valued Gaussian random variables of mean $ 0 $ and variance $ t $( cf. Random variable; Gaussian process). The convergence rate of the procedure is described by fluctuation results, large-deviation estimates, etc. (cf. also Probability of large deviations). The most elementary estimate concerns the standard deviation of the random error:

$$ { \frac{1}{\sqrt N } } \sqrt { {\mathsf E} \left | {u _ {0} ( W _ {t} ( x ) ) } \right | ^ {2} - \left | {u ( t,x ) } \right | ^ {2} } . $$

More generally, let $ L $ be an integro-differential operator (cf. also Integro-differential equation), and consider the evolution problem

$$ { \frac{d}{dt } } u ( t,x ) = Lu ( t,x ) , u ( 0,x ) = u _ {0} ( x ) . $$

From now on, it is assumed that there is a $ \mathbf R ^ {d} $- valued Markov process $ ( X _ {t} ) $ such that $ u ( t,x ) = {\mathsf E} u ( 0,X _ {t} ( x ) ) $( for a discussion on the assumptions that $ L $ must satisfy, see [a3], Chap. XIII, Sect. 3). Then, if it is possible to simulate independent copies of the process $ ( X _ {t} ( x ) ) $, one can approximate $ u ( t,x ) $ by

$$ { \frac{1}{N} } \sum _ {i = 1 } ^ { N } u _ {0} ( X ^ {i} _ {t} ( x ) ) . $$

For example, in neutronics, the process $ ( X _ {t} ) $ is the pair (position,velocity), and the velocity is a pure jump process (in this case, $ L $ is a transport operator) whose simulation is simple: one simulates the times at which the velocity changes and the new velocity after each change time independently (the change times have an exponential law).

Most often, the law of $ X _ {t} ( x ) $ cannot be simulated exactly, but it may be possible to approximate it suitably. For example, when $ ( X _ {t} ( x ) ) $ is the solution of a stochastic differential equation driven by a Brownian motion (in this case, $ L $ is a second-order differential operator), one can discretize the stochastic differential equation in time. The simplest discretization scheme (the Euler scheme) defines an approximate process $ ( X ^ {h} _ {t} ( x ) ) $ depending on the discretization step $ h $, whose simulation requires the simulation of independent Gaussian random variables only. The approximation formula is

$$ u ( t,x ) \simeq { \frac{1}{N} } \sum _ {i = 1 } ^ { N } u _ {0} ( X _ {t} ^ {h,i } ( x ) ) . $$

The global error of this procedure can be decomposed as the sum of

$$ u ( t,x ) - {\mathsf E} u _ {0} ( X ^ {h} _ {t} ( x ) ) $$

(the discretization error) and of

$$ {\mathsf E} u _ {0} ( X ^ {h} _ {t} ( x ) ) - { \frac{1}{N} } \sum _ {i = 1 } ^ { N } u _ {0} ( X _ {t} ^ {h,i } ( x ) ) $$

(the statistical error). Under appropriate assumptions, the discretization error can be expanded in terms of $ h $.

The numerical procedure can be adapted to the solution of problems in bounded domains. Dirichlet boundary conditions lead to the simulation of processes stopped at the boundary. The simulation of reflected diffusion processes (cf. Diffusion process) solves parabolic or elliptic partial differential equations with suitable Neumann boundary conditions.

Let $ L ^ {*} $ be the formal adjoint of $ L $. Note that the probability distribution $ \mu _ {t} ^ {x} $ of $ X _ {t} ( x ) $ is a solution to

$$ { \frac{d}{dt } } \mu _ {t} ^ {x} = L ^ {*} \mu _ {t} ^ {x} $$

in the sense of distributions (cf. Generalized function). The solution $ \mu _ {t} ^ {x} $ can be approximated by the empirical distribution

$$ { \frac{1}{N} } \sum _ {i = 1 } ^ { N } \delta _ {X _ {t} ^ {i} ( x ) } . $$

This procedure is used in neutronics.

To certain non-linear partial differential equations one can associate processes $ ( X _ {t} ) $ which are not Markovian, but are such that the pairs $ ( X _ {t} , \mu _ {t} ) $ are Markovian, where $ \mu _ {t} $ is the probability law of the random variable $ X _ {t} $. Besides, in some cases, the law $ P ^ {X} $ of the process $ ( X _ {t} ) $ can be described by means of the theory of propagation of chaos (see [a4] and the references therein). That means that the law $ P ^ {X} $ is the limit in probability, when $ N $ goes to infinity, of the empirical distributions of the paths of $ N $ particles in a mean field interaction. The interaction between the particles is related to the non-linear differential operator of the partial differential equation under consideration. If one can simulate the stochastic particle system, then an algorithm can be built from the formula

$$ \mu _ {t} \simeq { \frac{1}{N} } \sum _ {i = 1 } ^ { N } \delta _ {X _ {t} ^ {i} } , $$

where $ \{ X ^ {i} _ {t} \} $ are the locations at time $ t $ of the interacting particles. For example, consider an interaction kernel $ K $, a family of $ N $ independent Brownian motions, and the system

$$ X _ {t} ^ {i} = X _ {0} ^ {i} + \int\limits _ { 0 } ^ { t } { { \frac{1}{N} } \sum _ {j = 1 } ^ { N } K ( X _ {s} ^ {i} - X _ {s} ^ {j} ) } {ds } + \nu W _ {t} ^ {i} . $$

The time discretization of such systems permits one to approximate equations such as Vlasov–Poisson–Fokker–Planck equations (cf. Vlasov–Poisson–Fokker–Planck system), the viscous Burgers equation (cf. Diffusion equation; Turbulence, mathematical problems in), the incompressible Navier–Stokes equation for the vorticity in $ \mathbf R ^ {2} $( cf. Navier–Stokes equations), etc. In the latter case, $ K $ is the Biot–Savart kernel, and the algorithm is the Chorin random vortex method (see [a2] for numerical considerations). Similarly, the simulation of interacting branching processes (cf. Branching process) permits one to approximate the solutions of convection-reaction-diffusion equations. One difficulty of the mathematical study of these methods is to show that the statistical errors (on the distribution function of $ \mu _ {t} $, e.g.) have standard deviations of order $ {1 / {\sqrt N } } $ although the particles are dependent. For theoretical aspects related to probabilistic algorithms for partial differential equations, including the Boltzmann equation, see, for instance, [a1] and the references therein.

References

[a1] C. Graham, T. Kurtz, S. Méléard, P. Protter, M. Pulvirenti, D. Talay, "Probabilistic models for nonlinear partial differential equations" , Lecture Notes in Mathematics , 1627 , Springer (1996)
[a2] "Vortex methods and vortex motions" K.E. Gustavson (ed.) J.A. Sethian (ed.) , SIAM (1991)
[a3] J. Jacod, "Calcul stochastique et problèmes de martingales" , Lecture Notes in Mathematics , 714 , Springer (1979)
[a4] A.S. Sznitman, "Topics in propagation of chaos" P.L. Hennequin (ed.) , Ecole d'Eté de Probabilités de Saint-Flour XI (1989) , Lecture Notes in Mathematics , 1464 , Springer (1991) pp. 165–251
How to Cite This Entry:
Monte-Carlo methods for partial differential equations. Encyclopedia of Mathematics. URL: http://encyclopediaofmath.org/index.php?title=Monte-Carlo_methods_for_partial_differential_equations&oldid=47896
This article was adapted from an original article by D. Talay (originator), which appeared in Encyclopedia of Mathematics - ISBN 1402006098. See original article