Python as a powerful calculator

We start our tour of Python with its built-in data structures. The simplest data type we can think of are numbers. Python has two data types to deal with them: int and float. The type int is used to represent integer numbers, while float represents decimal point numbers. We could as an example check the types of 3 and 3.0.

print("type of 3: ", type(3))
print("type of 3.0: ", type(3.0))


Exercise 1. Use shift + enter to run the code above in the box below.



In [ ]:



The function type() you have just used is a very important construct that highlights the fact that everything in Python is an object of a particular type. The type of an object determines what operations we can perform on them. Of course, our purpose in programming is to manipulate data. For numbers in particular, we expect to be able to compute all the familiar algebraic operations we learned at school: sum, subtraction, multiplication, division and exponentiation. The main algebraic operations are implemented in Python in the following way:

• Sum: +
• Subtraction: -
• Multiplication: *
• Division: /
• Exponentiation: **
• Integer division: //
• Modular algebra: %

In the interactive modes, be it a Jupyter notebook or an IPython console in Spyder, Python can indeed be used as a powerful calculator.

Exercise 2. Compute the expression $\sqrt[4]{\frac{765\cdot(587 + 13987)}{4^{15}}}$.



In [ ]:



A few relevant points are worth to note when computing expressions like the one above. We all carry from school the notion of operator precedence. Thus, products should be solved before sums, and exponentials before everything else. These common sense rules remain true in Python. Evidently, we can always change the order of precedence with parantheses, just as we do in real life. You can check in the following link the complete list of operator precedence in Python.

Another important point to note is that some operations output data of a different type from their inputs, as we see below.

print("type of 1: ", type(1))
print("type of 1/1: ", type(1/1))


Exercise 3. Run the last commands in the box below.



In [ ]:



In the particular case of the division operator, this has not always been the case. In Python 2, the division operator applied to two integers yielded their integer division as result. Integer division is accomplished in Python 3 with the operator //. The modulo operator, which yields the remainder of integer division, is %.

Exercise 4. Compute 1//2 and 8%3.



In [ ]:



Variables and assigment

When programming, we will rarely perform these algebraic operations directly on numbers like we do on a calculator. We use variables instead. In Python, the assignment operator is =. The expression a = 3 means that we are assigning the value 3 to the variable a. The piece of code below shows how things are usually done.

a = 3
b = 6
product = a*b
print(product)


Exercise 5. Run the code written above.



In [ ]:



A brief comment is due on assignments in Python. We should think of the command a = 3 as binding the name a to the int object 3. Names in Python live in namespaces, which can be thought of as maps from names to objects. The details of namespaces and scopes will become more clear as we increment our use of the language. At this point, specially when dealing with numerical types, it won't hurt to think of the object and its name interchangeably. However, you should always have in mind that a variable is just a name and nothing more. Ultimately names and the things they refer to are distinct entities, as is known since the days of Juliet. Nevertheless, you should be judicious when choosing variable names. Good names are usually descriptive and instantly inform someone who is reading the code what is the purpose of the corresponding variable. In this series of notebooks we try to abide to this common sense principle.

As a final comment on this last piece of code, note that we have written the statement print(product) to print the result of the product to the screen. In the interactive mode, we might simply have typed a*b after assigning a and b since the results of operations are automatically printed. In the script mode though, we have to use the print statement.

Functions and plotting

We now introduce two more concepts that will be recurring during this course: those of function and plotting. Frequently we encounter situations in which we know the relation between inputs and outputs, but we don't know precisely what the input is. The notion of function then comes handy. Suppose for example that we want to evaluate the function $f(x) = 3x^2 + 2x - 4$ for several values of $x$. We can easily define this function in Python:

def f(x):
return 3*x**2 + 2*x - 4


After the command above is run, you can use expressions like y = f(x) to assign the result of the expression 3*x**2 + 2*x - 4 to the variable y.

Try this below.

Exercise 6. Define the function $f(x) = 3x^2 + 2x - 4$.

Use your code to compute $f(4298)$. Assign the result to the variable y and print y to the screen.



In [ ]:



The piece of code above illustrates how functions are implemented in Python. The general structure

def f(x):
indented code
return something #this is optional!


is always used to define a function $f$ that depends on a parameter $x$. The (optional) command return causes the function to exit as soon as it is reached and passes the subsequent expression to whatever called the function. One important point to note is the colon : and the indentation. Everything we write inside the definition of the function and up to the return statement has to be indented. These indentation rules will be recurring whenever we change scope, which in Python is indicated by the colon. The concept of scope is intimately related with namespaces. It is yet another concept that is best learned by practice and will surely become clearer as we advance. As to the actual number of spaces used for indentation, it is yours to choose. You must only stick to your choice (i.e., always use the same number of spaces) when changing scope through your code. For readability we suggest indenting four spaces.

Finally, we can plot the function $f$ that we have just defined. To do this we will first import the pyplot module from the matplotlib library. This is achieved through the following line of code.

import matplotlib.pyplot as plt


The statement above finds, loads and initializes the module pyplot that exists in the matplotlib library. In particular, we have specified that the imported module shall be referred to as plt, so that we don't have to type the full expression matplotlib.pyplot whenever we want to use it. Thus we generally do something like

import long_module_name as short_name


so that we don't have to write the full module name whenever we use it. More details can be found in the references on importing. We are ready to do our plot.

import matplotlib.pyplot as plt
X = [(-10 + i*0.1) for i in range(200)]
Y = [f(x) for x in X]
plt.plot(X,Y, color = "r")
plt.xlabel("x")
plt.ylabel("f(x)")
plt.title("My first plot!")
plt.show()


The last lines have introduced many concepts that we shall be using again and again. We will not go in detail at this point, but here goes a brief description of what is going on: both X and Y (capitalized! Python is case sensitive) are examples of iterables, which are data structures similar to lists. In other words, they have some internal structure. You can check their types as we have done before. Note in particular how the lists X and Y were created. This very practical method of creating lists is called list comprehension. The subsequent lines use some of the functions provided by the matplotlib.pyplot module. The lists X and Y are both a discrete collection of numbers. You can look at their contents using the print() command. We chose their length to be 200 so that the resulting plot looked smooth, but you should always remember that there is no such thing as a real number in a computer! The plotted curve is actually polygonal. The code specifies that the plot should be red ("r") and also includes lines that provide labels to the x and y axes of our graph, as well as its title. The last line, plt.show(), is necessary to visualize our plot. In the interactive mode, it should generate a picture with the plot as output. In the script mode, a new window with the plot should pop up. You should try to put all these lines of code together in Spyder to see what happens in the script mode.

Exercise 7. Run the code written above to plot $f$.



In [ ]:



Exercise 8. The Lennard-Jones potential.

The Lennard-Jones potential is the simplest way to model interactions between electrically neutral molecules. It is usually written as

$$\frac{V_{LJ}(r)}{\epsilon} = \Big[\Big(\frac{r_m}{r}\Big)^{12} - 2\Big(\frac{r_m}{r}\Big)^6\Big]$$

where $\epsilon$ is the depth of the potential well and $r_m$ is the distance where the potential is a minimum. Note that both sides of this expression are dimensionless. Define $x = \frac{r}{r_m}$ and $y = \frac{V_{LJ}(r)}{\epsilon}$ and then plot the graphic of $y$ versus $x$.



In [ ]:



The math library

In the next exercises we will play a little bit with the math library. It contains mathematical functions and constants commonly found in science and engineering such as $\sin(x)$, $\exp(x)$ as well as constants like $\pi$. Although we will usually prefer numpy for our numerical tasks, it is a good idea to be comfortable with the math module too.

Exercise 9. Compute the area and perimeter of a circle.

Below we start the definition of two functions, circle_area(radius) and circle_circumference(radius). Use the math library implementantion of $\pi$.



In [6]:

import math
'''
Returns the area of a circle
'''
###Substitute the pass command with the answer
pass

'''
Returns the circumference of a circle
'''
###Substitute the pass command with the answer
pass



Exercise 10. The Yukawa potential.

The Yukawa potential is a generalization of the Coulomb potential that is exponentially damped by a mass term. It can be scaled to be written dimensionlessly as

$$\frac{V_{Yukawa}(r)}{g^2km} = - \frac{e^{-kmr}}{kmr}$$

where $m$ is the mass of the particle and $k,g$ are scaling constants. Define $x = kmr$, $y =\frac{V_{Yukawa}(r)}{g^2km}$ and then plot the graph of $y$ versus $x$. You can also check how the Yukawa potential changes when you fix $k,g$ and vary $m$. See exercise 13 below and this picture in Wikipedia.

Exercise 11. Real roots of quadratic polynomials.

In this exercise we implement in Python a solution to the quadratic equation,

$$ax^2 + bx + c = 0.$$

We shall assume that the equation has real solutions, which means that its coefficients are such that

$$b^2 -4ac \geq 0.$$

Fill the code below. You might want to use the square root function in the math module. After you have finished, use your code to find the golden ratio as a solution to the equation $x^2 - x - 1 = 0$.



In [17]:

'''
Returns the roots of the quadratic equation ax**2 + bx + c = 0. Assumes b**2 - 4ac > 0.

a: Coefficient multiplying x**2
b: coefficient multiplying x
c: independent term
'''
pass



Exercise 12. Heron's formula.

A classical result from euclidian geometry is Heron's formula, which computes the area of any triangle once its sides are known. It states that the area $A$ of a triangle with sides $a$, $b$ and $c$ is given as

$$A = \sqrt{s(s-a)(s-b)(s-c)}$$

where s = $\frac{a+b+c}{2}$ is the semiperimeter. Implement Heron's formula below.



In [12]:

def semi_perimeter(a,b,c):
"""
Returns the semi_perimeter of a triangle.

a: Side of the triangle.
b: Side of the triangle.
c: Side of the triangle.
"""
pass

def heron(a,b,c):
"""
Returns the area of a triangle.

a: Side of the triangle.
b: Side of the triangle.
c: Side of the triangle.
"""
pass



Exercise 13. Shannon entropy.

In information theory, Shannon entropy is a measure of the amount of information generated by a stochastic source of data. Given a discrete random variable $X$ with probability mass function $P(X)$, its corresponding Shannon entropy will be

$$H(X) = E[-\ln P(x)] = -\sum p(x)\ln p(x)$$

Define and plot the function $f(p) = p \ln p$.



In [ ]:



Exercise 14. Plot several cosine functions with different angular frequencies.

To make our task a little easier, we will exploit Python's object orientation. Since everything in Python is an object, we can move things around without much difficulty. In this exercise, we will define a function cosine(omega) that returns a cosine function with given angular frequency.

import math
def cosine(omega):
"""
Returns the function cos(omega*t).

omega: The angular frequency.
"""

def cosine_omega(t):
"""
Returns the number cos(omega*t).

t: A float number.
"""

return math.cos(2*math.pi*omega*t)

return cosine_omega


The idea behind the function cosine(omega) is to implement the family of functions $f_\omega(t) = \cos(\omega t)$ that are labelled by the parameter $\omega$. Thus we map $f_\omega \mapsto$ cosine(omega). Run this code in the cell below.



In [14]:



Now you can use the function cosine(omega) defined above to plot the functions $\cos_k(t) = \cos(2k\pi t)$ for $k = 1,3,5,7$ and $t$ in $[0,1]$. Use different colors and different markers for each function. Label your axes, put legends identifying each function and give your plot a title. Finally, save your plot as "cosines plot.pdf". To make your life easier, the function $\cos_i$ is already worked out. This documentation on pyplot might be helpful.



In [15]:

import math
import matplotlib.pyplot as plt

cos_1 = cosine(1)

### Uncomment and complete the following

#cos_3 =
#cos_5 =
#cos_7 =

X = [i*0.01 for i in range(100)]
Y1 = [cos_1(x) for x in X]

plt.plot(X,Y1, color = 'red', linestyle = 'solid', label = '$k = 1$')
plt.legend()
plt.show()






Exercise 15. Plot several normal distributions.

Your plot should show several normal distributions $\mathcal{N}(\mu,\sigma)$ with same mean $\mu$ and different standard deviations $\sigma$. Take for example $\mu = 0$ and $\sigma = 1, 2, 4$. Follow the same lines as the previous exercise.



In [ ]:



Exercise 16. Subplots of Legendre Polynomials.

The Legendre polynomials often appear in physics, in particular related to solutions of Laplace's equation. In the wikipedia link given you will find the first eleven of them, as defined in the $[-1,1]$ closed interval. Try to reproduce the picture in wikipedia exactly. After this, make subplots of each polynomial separately.



In [12]:



Exercise 17. Scatter plot of pi.

Below you are given a list that contains the first 101 digits of $\pi$ taken from here. Imagine that $\pi$ is actually a map $f:\mathbb{N} \rightarrow \mathbb{N}$, such that $f(0) = 3$, $f(1) = 1$, $f(2) = 4$ and so on. Make a scatter plot of this map.



In [46]:

pi = list('31415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421170679')



Exercise 18. Histogram of $\pi$ digits.

Plot a histogram of $\pi$ digits.



In [ ]: