GEOL 5470

Homework/Labwork 1                                                                                                     Humphrey

This is mainly to get you working with Python and Spyder and using some of the processing and visualization techniques.  For anything that doesn’t make sense, probably the fastest way to get information is to Google.  For example to see how a python for loop works, just google ‘python for loop’.  The 2 main web resources are ‘docs.python.org’ for basic tutorials and documentation, and ‘stackoverflow.com’ for answers to weird behavior that you don’t understand.  There are LOTS (or maybe too many) of resources out there. The first few steps below are actually covered in the PYTHON primer introduced in previous notes. There are several example programs on our web page that do this homework, but I would like you to try to do it without just copying. You will only learn by making mistakes.

1.       Open Spyder, and make a ‘new’ file in the editor. Add the lines ‘import numpy as np’ and ‘import matplotlib.pyplot as plt’, this will allow you to use the math of numpy methods as ‘np.method’ and the plotting abilities of matplotlib as ‘plt.method’ in your program.

2.       First project. Create a vector which is a sine wave. Use the numpy method np.linspace(start,stop,number of elements) to make a 1D array of numbers (call it x, we will use it again later), and then use the np.sin() method to create a list of sine values, assigned to the variable w. Scale start and stop so that w spans several oscillations of the sine wave.  Numpy trig functions use radians, not degrees ( radian = degrees * pi/180). The trig functions accept lists or arrays, and return arrays. Comment your program with the # character!  Run the program, and check in the Ipython Console that it ran, and that w looks like a vector of numbers.

Note that these 2 lines are an example of vectorization or, using the pythonic term, broadcasting. The first line creates a long list of number (actually a 1D array), while the second line that invokes the sin() method operates on the entire array. In general, numpy methods always operate on all the elements of any array you give them.  Numpy arrays print out in the console with ‘[‘ type brackets around each row.

3.      Second project.  Plot the Sine of x.  You can just type: plt.plot(x,w) to get a bare-bones plot. (add plt.show() to actually see your plot).  Just to have some consistency, and to allow us to modify our plots in future, we will be a little more verbose than necessary.  We will first make a figure window on which to plot: fig = plt.figure().  Then make a set of axes: ax = fig.add_subplot(1,1,1).  Note that the first line created a figure object, and in the second line we exploit the fact that figures have methods, and one of those methods is called ‘add_subplot’, which makes one or more plot axes on the figure.  Lookup matplotlib add_subplot if you want to know more (you should find the matplotlib documentation).   Now to actually plot your data on the axes: ax.plot(x,w), which produces a minimal plot.  However, although the plot has been made, you need to tell python what to do with the plot, so add a last line: plt.show(). (Typically, the ANACONDA version automatically defaults to including plt.show(), so your plot will be visible without it, but some other programs don’t, so it is good habit)

Make sure that ‘tools/preferences/ipython console/graphics’ is changed from ‘inline’ to ‘automatic’ to get the plot to show in a separate window.  The plots sometimes end up BEHIND the Spyder GUI! On the desktop, so look for it.

4.  (A little advanced) Modify your plot by adding an x axis in degrees, include title and x, y labels. Add line color, line style and markers.  This can be done 2 ways.  You can directly edit the plot, using the tools at the top of the plot.  This is good for quick and dirty changes.  The best method however is to use the ‘set’ methods, such as: ax.set_title(“Plot of the Sine function”).  Modifying the plot in your program has to occur after you have made the axes.  Options for plt.set can be found by searching for matplotlib axes class.  As you can see, you can do just about anything!

Matplotlib is a HUGE set of methods that can make very good 2D and pseudo 3D plots.  It is a steep learning curve, and also a big step into Python.  But plotting is core to using numerics for science so I think the jump is worthwhile.

Part B

Now to start to learn the single most important part of most numerical methods: to get the computer to do math for you.  This will expand our abilities in Python.

1.       Add a cosine curve to your plot, in a different color, using the x same numbers as the sine curve. Note, cosine is the first derivative of sine, so you have plotted a function and its derivative. But now will try to calculate the derivative numerically.

2.       Make a vector of the differences between consecutive numbers in sin(x), and divide each difference by the difference in x (delx) between each value (this should be a constant).  (Note you can ask numpy to do this differencing by calling np.diff(w), Most simple tasks have already built methods in numpy, but for now we are going to do it the hard way).  To take this difference in python, you will need to write a ‘for’ loop.  The syntax of a ‘for’ loop is: a line containing ‘for x in y:’ followed by any number of python statements, all consistently indented. Note the ‘:’ at the end of the line, important!  There is no ‘end’ statement, only a reduction in indentation.  The ‘x’ stands as a dummy variable for each value of the variables in y.  Unlike many languages, y can be anything that can be indexed or iterated.  Typically y is a list of numbers, but it can also be a string, or a list of colors etc.  y is often a generator of a list, such as np.arange(start,stop,step), which is particularly useful to take integer steps and to mimic ‘for’ loops in languages like C and Fortran.  So our ‘for’ loop to create a ‘diff’ vector:

# assuming you have an array of numbers ‘w’ 

diff = np.zeros(len(w)-1) # create 1D array as a place to put our calculated differences

for i in np.arange(0,len(w)-1): # loop thru the numbers and calculate the numerical difference between each

    diff[i] = w[i+1] – w[i]

diff = diff/(x[1]-x[0]) # turn the difference array into a ‘numerical derivative’ by dividing by the radians between numbers

This is our first ‘real’ python numeric code, so it needs some explaining.  Note, ‘copy and paste’ from HTML (the language of the web page) will often include non-printing HTML characters that make python complain.  For now it is best to type these lines into your program (in the above code the ‘-‘ sign is not an ASCII character in HTML.)

First line uses np.zeros(number of points) to create an array to put our result.  It is one shorter than the length of w. The second line is the start of the for loop, which uses the ‘generator’ np.arange() to make a list of integers from start to 1 less than stop.  The generated list is used as an index ‘i’.  The index ‘i’ can be used to select a single value from the array w.  There is only one statement in the loop, which takes the w value at i+1 and subtracts the value at i, putting the result in ‘diff’ at the I location.  At the end of the loop, reduce the indent to indicate the loop end, and we use a vectorization to divide each of the diff values by the uniform step size in x.

3.       Plot the resulting difference vector, on the same axes as all the above, with the points plotted at ½ way in between the w points (the difference vector is 1 point shorter than the w vector, so it will all fit).  You may run into problems because you will want to plot:  plt.plot(x,diff), but x and diff have different lengths.  The solution is to plot: plt.plot(x[:-1],diff).

Note on a fairly extensive topic:  the term x[:-1] above is a python basic method that works on lists of things, including arrays.  It is referred to as indexing and/or slicing an array.  At its simplest, w[i] just selects a single value, note the first value has index=0.  There are some useful special cases:  w[-1] selects the last value, w[-3] selects the 3 from last etc.  When selecting more than one value, this is called slicing.  In slicing, the symbol ‘:’ means ‘all the values inbetween’.  So w[:] is the same as just w (everything).  And w[0:-1] means include from the start (0), all inbetween, but stop at the last value (slicing includes the start but not the stop).  This can be shortened to w[:-1] since the (0) is the default start.  The default stop is the end, so w[5:] is a slice from 5 to the end.  You can go crazy with slicing, using negative steps etc, (see the docs) but you can also make your code hard to read if you get too fancy.

4.       Compare the difference vector to the cosine curve.  With any luck they should be very similar: in other words we have developed a technique for differentiating a function.  This works for any function. 

5.       To do a little numerical experimenting, try changing the step size in x, you will see that when you take big steps the quality of your differentiation will decrease.

6. If you find this easy, then try putting a legend on the plot.

7.       Try differentiating a polynomial, see if x^3 differentiates into x^2 /3.