# Preface

A numerical approximation to the solution of a linear elliptic partial differential equation at a single point. Simulation of the motion of a random particle may be used to approximate the solution to linear elliptic equations.

Return to computing page for the first course APMA0330

Return to computing page for the second course APMA0340

Return to Mathematica tutorial for the first course APMA0330

Return to Mathematica tutorial for the second course APMA0340

Return to the main page for the first course APMA0330

Return to the main page for the second course APMA0340

Return to Part VI of the course APMA0340

Introduction to Linear Algebra with *Mathematica*

## Glossary

# Monte Carlo for Elliptic Equations

Simulation of the motion of a random particle may be used to approximate the solution to linear elliptic equation. In particular, consider the two-dimensional Dirichlet boundary value problem:

*L[ ]*defined by

*A, B, C, D, E, F*are all functions of {x,y,t}. Here

*R*is two dimensional domain on ℝ² with smooth boundary ∂

*R*.

The steps for this method are straightforward. First, approximate the given elliptic partial differential equation by a finite difference method. Rewrite the finite difference formula as a recursive function for the value of the unknown at any given point. Then interpret this recursive formula as a set of transition probabilities that determine the motion of a random particle.

Now, write a computer program that will allow many (say *K*) particles
to wander randomly around the domain of interest, based on the
transition probabilities found from the difference formula. Simulate
particles one at a time, with every particle starting off at the same
point (say the point $\bfz$).

If the boundary data are not necessarily of the Dirichlet type (i.e., Neumann or mixed boundary conditions) then, when the particles reach the boundary, they will be given a finite probability to leave the boundary, and re-enter the domain of the problem. If the particle leaves the boundary, then continue the iteration process. If it does not leave the boundary, then the value at the boundary is stored away, and a new particle is started off at the point~$\bfz$.

The simulation is finished after all $K$ particles have been absorbed into the boundary. If the original elliptic equation was homogeneous, then an approximation to the solution, at the point $\bfz$, will be given by the average of all the values obtained (recall that when the particles stop at the boundary they obtain a value).

step = 0.10; SUM = 0.0; disp("Number of particles Approximation") for i = 1:10000 X = [2;0]; r = norm(X); while r <= 3 && r > 1 X = X + step*sign(rand(2,1) - .5); r = norm(X); end if r <= 1 % When a particle hits the boundary, sum the value SUM = SUM + 4; end if r >= 3 SUM = SUM + 6; end if mod(i,1000) == 0 fprintf("%10d %21.3f \n", i, SUM/i) end end

Return to Mathematica page

Return to the main page (APMA0340)

Return to the Part 1 Matrix Algebra

Return to the Part 2 Linear Systems of Ordinary Differential Equations

Return to the Part 3 Non-linear Systems of Ordinary Differential Equations

Return to the Part 4 Numerical Methods

Return to the Part 5 Fourier Series

Return to the Part 6 Partial Differential Equations