Hello, dear friend, you can consult us at any time if you have any questions, add WeChat: THEend8_
Prerequisites:
MATH1715, MATH1725, some programming experience in a computer language like Python or R.
Preamble
During this project you will get a taste of doing original mathematical research. This involves being proactive in planning your project and ?nding ways to make things work. You will work partly in a team but mostly independently, to produce a group presentation and an individual report on your project.
We will begin by studying some topics together, and discussing some introductory homework tasks to gain experience of the topics. Later (for perhaps the last 60% of the project), each person will choose a different topic from the list in section 5 to investigate in greater depth, or may propose another topic to their supervisor.
All projects will require some investigation of the background theory, and some computational investigation of speci?c problems, but the balance between these is up to you. Your report will present all of your investigations, both the introductory ones and the more substantial one, explaining the background, results and conclusions in a style similar to that of a textbook, so that they could be understood by a maths undergraduate not familiar with the project.
1 Introduction to Monte Carlo methods
Most computer programming languages have a command which allows the user to generate random, or pseudorandom, numbers. For example, the “Python” language has a set of standard commands for generating random numbers (either integers or real numbers) from a wide range of standard distributions, such as the uniform and normal distributions [1]. The aim of this project is to investigate computational methods, known as Monte Carlo methods, which can solve a variety of mathematical and physical problems, to arbitrary accuracy, by using random numbers.
Random numbers can be used for statistical sampling (also called “Monte Carlo sampling”). A common application of this is in numerical integration, in which the function to be integrated is sampled at random points within the domain of integra- tion, in order to estimate the value of the integral [2, 3]. For integration in one dimen- sion, other numerical integration methods (such as Simpson’s rule) are usually more ef?cient than Monte Carlo Methods. But for multidimensional integrals, Monte Carlo methods become competitive, and for very large numbers of dimensions these are the only practical method for numerical integration.
Homework task 1:
The following are to be completed before the second supervised meeting. Please create presentable results that you can share on your screen (or upload to the TEAMs channel) for the second meeting.
You may choose to work in any standard programming language to which you have access, such as Python, R, C or Java. During the supervised sessions we will study examples in Python.
? Preparation. You should begin by making sure that you can edit a program ?le and execute the program (after ?rst compiling it, if you are using a compiled lan- guage). Information on how to obtain Python has been sent by email. The best solution is to install a suitable programming language on your own computer. For Python we recommend installing the anaconda package which includes the integrated development environment (IDE) “Spyder” . “R Studio” provides a sim- ilar functionality for R if you prefer to use that.
? When learning a programming language or using an IDE for the ?rst time, it is traditional to begin by writing a program to output the text “Hello World” . Make sure that you can do so. (This task is too simple to be included in the ?nal report.)
? Run the Python code provided in the session (or your own version of it in your chosen programming language) to evaluate the mean and variance of a random sample drawn from the probability distribution(s) discussed in class. Then adapt the code to investigate different parameters and/or distributions as discussed. Check that the results make sense to you. If they do not, it may indicate a bug in the program. Edit the program until you understand what is happening.
? For your report and presentation you will need to be able to present your re- sults in the form of graphs and ?gures. In Python there is an extensive library of graphing functions called Matplotlib. One of the example Python codes uses this to plot a histogram of sampling from a uniform distribution, but this pack- age can do many other types of plots (for more examples and tutorials go to https://matplotlib.org). Create a ?gure or ?gures to show how the number of samples affects the accuracy of capturing a distribution.
2 Monte Carlo Integration
A good reference for this section is chapter 7.7 of Ref. [4], which is available online at http://numerical.recipes/book/book.html or there are print copies in the library.
A common task in mathematics is the evaluation of integrals. Although in the problems set in undergraduate mathematics courses integrals can be performed an- alytically, this is rarely true in applications. Consequently we need numerical methods to approximate integrals as sums, so that
where we take a weighted sum of the values of the function at N points xi. Intradi- tional quadrature schemes such as the trapezium and Simpson’s rules [4] the posi- tions of the points are predetermined, e.g. equally spaced points.
However, in Monte Carlo integration the points are chosen randomly from a distri- bution. Since the mean value off can bede?ned as
we can determine the value of the integral by sampling f to ?nd its mean value. Thus if we draw N points, x1 , . . . xN from a uniform distribution in the interval [a, b] then we can estimate the integral as,
Furthermore we can also estimate the error, by taking the square root of the mean squared error as
where σ 2 is the variance of f over the interval. This means that the accuracy of the estimate improves as the square root of the number of points sampled.
More generally if points are chosen from a probability distribution p(x) then
Homework task 2: For a general function f (x) obtain a formula for σ and hence an error estimate. How does this compare with quadrature schemes such as the trapezium rule or Simpson’s rule? Write a program to estimate R 0 2 cos2(x) dx using Monte-Carlo sampling. Compare the answer with a classical method such as the trapezium and Simpson’s rule. Is Monte-Carlo sampling a good way to estimate this integral?
2.1 Multi-dimensional Integration
The above example shows that there are better ways of integrating in one-dimension than Monte-Carlo (MC) methods. However, we shall discover that in multi-dimensional integration MC methods have their advantages.
Suppose we now consider the integration of the function f (x) over a D-dimensional domain x ∈ V. As shown in the vector calculus course this can be computed via nested one-dimensional integrals as
where V is contained in the hypercuboid (a, b) and we de?ne f to be zero for points outside V. If we use a quadrature scheme with n points for each of the integrals we require a total of N = nD function evaluations. Now unless V ?ts precisely to the boundaries of the domain, f will be discontinuous. This means the error in each quadrature integral will scale as 1/n (since the higher order convergence of the trapezium and Simpson’s rules depend on function having continuous derivatives). Thus the total error will scale as N —1/D. So when D is large the convergence is very slow. Even where we have a rectangular domain over which the function is differ- entiable the convergence rate we decrease as the number of dimensions increases. This is the “curse of dimensionality” .
Now suppose we use Monte-Carlo sampling. The integral can be estimated as
Hence the error scales as N—1/2 irrespective of the number of dimensions.
Homework task 3:
As an example let us consider calculating the area of the unit circle by integrating the function
over the domain V = {(x1, x2) : |x1| < 1, |x2| < 1}. The answer is of course π . Write and execute your own program to estimate the integral using Monte-Carlo sampling. By plotting a graph, examine how the error changes with N. Run repeated calcula- tions for particular values of N and compare the variance with the error estimator.
How many samples would you need to calculate π to 3, 5 and 10 decimal places?
As usual, present your program and results at the next supervised meeting.
Extension task i:
Modify your programme to calculate the volume of a unit sphere and higher dimen- sional hyperspheres. The analytical formulas for the answers are given at http://thinkzone.wlonk. com/MathFun/Dimens.htm. How does increasing the dimension modify the convergence?
Extension task ii:
Examine what happens if you replace f with the function
where In particular look at what happens for the cases —2 < Q < 0 when the integrand is singular at (0, 0) but integral remains convergent. Also see what happens at Q = —2 when the integral has a logarithmic divergence.