Phys690: Homework Assignment 1

==> Due on Monday 2/19 in class <==

(1). Error analysis and the central limit theorem.

File data contains a set of results from the pi calculation we discussed in class. It has n=2500 numbers, each a measurement of pi with m=10000 points thrown in the square. Here we will analyze this data using Maple.

Download the data file, e.g., as 'data' in a directory that you make for this homework. On the linux machines in Physics (see 'things to note when computing' and 'Maple quick reference page'), you can start a regular session by typing maple or the X-window version by xmaple. (With the latter, display must be set properly.)

Attached is a list of Maple commands which will do most of the analysis. (They can also be found in a separate file commands, which you can download.) Please feel free to just run these to get going, but please be sure that by the end of this exercise you understand what the code is doing.

1. What is the average of these measurements? What is their standard deviation?
2. Estimate the statistical error of the average above. How does our final result compare with the exact value of pi?
3. If, instead of the entire set, we take a fraction, say, the first 100 numbers, and calculate their standard deviation, how will it compare with the standard deviation computed in 1? If you are not sure, do the calculation and/or study the histograms. How would the statistical error of the average of these 100 numbers compare with that of the entire set?
4. As you see, it takes a little time for Maple to read in a list of 2500 elements. In order to reduce this problem, let us imagine a perhaps more civilized approach to generate the data file: Rather than recording every number, we average every 50 numbers and only write out the average. This way the data file would have only 50 entries, with the first one being the average of the first 50 numbers in our current file, and the second one being the average of the next 50 numbers, and so on. With this new file, what would be the answers to 1? What about 2? Justify your results with either mathematical arguments or numerical experiments (or better yet, both). If you choose to do numerical experiments, you can either write a program to actually produce the data file as described, or manipulate inside Maple to generate the 50 numbers and then directly use them to get the answers.

Maple commands: ('#' indicates comments.)

with(stats):             # These 2 lines read in Maple packages for histogram; 
with(stats[statplots]):  # ':' at the end of a statement suppresses output;
                         # ';' allows it.
aa := readdata(data):    # Read data into a list called 'aa'.
ndim := nops(aa);        # number of elements in list 'aa'.
histogram(aa,color=yellow); 
                         # Just to see what the distribution looks like;
                         # note graph size can be adjusted.
xavg := sum(aa[i],i=1..ndim)/ndim;
xsigma := sqrt(sum(aa[i]^2,i=1..ndim)/ndim - xavg^2);
error := xsigma/sqrt(ndim);
#
bb := aa[1..100]:        # Assigns to list 'bb' the portion of 'aa' 
                         # between the 1st and the 100th elements.

(2). Drawing a Fern.

As we have learned, natural processes that are subject to chance influences can produce objects of high regularity. To demonstrate this, write a program to iteratively generate points in two-dimensional space using the following rules:

    (x_{n+1},y_{n+1}) = 
                    
  / (0.5,0.27*y_n),                                                  with 2% probability;
  |                
  | (-0.139*x_n + 0.263*y_n + 0.57, 0.246*x_n + 0.224*y_n - 0.036),  with 15% probability;
 <
  | (0.17*x_n - 0.215*y_n + 0.408, 0.222*x_n + 0.176*y_n + 0.0893),  with 13% probability;
  |
  \ (0.781*x_n + 0.034*y_n + 0.1075, -0.032*x_n + 0.739*y_n + 0.27), with 70% probability.

Start from an initial point (x_1,y_1)=(0.5,0.0). Carry out the iteration at least 30,000 times and plot all the data you obtain (as points) in an x-y plot.