Anne Kværnø
Date: Aug 18, 2022
A Jupyter notebook is composed of cells. These can contain either text (like this one) or code (as demonstrated below). As soon as you have copied the Jupyter note file to your own domain, you are free to use it and change it as you like. You are encouraged to experiment with the codes. You can not make any harm, if your note is completely messed up, just make a new copy.
If you have installed the Anaconda Python distribution as suggested above,
you can use Anaconda Navigator
to start the Jupyter notebook application. If you prefer the command line
in a terminal, you can just type jupyter notebook
.
There is plenty of documentation regarding the Jupyter Notebook which is part of the Jupyter project, see for instance
Also, the JupyterLab is a modernized alternative to the Jupyter Notebook. It can be launched from the Anaconda Navigator as well.
As soon as your first notebook is open, we recommend that you to make yourself more familiar with the Jupyter environment. Take a look at the tool menu, and see what is available. There are also many keyboard shortcuts (see Help -> Keyboard Shortcuts); learning some of those will make your work faster
Let us start: Below, you will find a short example of code. You can run this either by pressing either one of the following short cuts:
shift+enter: executes cell and takes you to the next cell (if there is one)
ctrl+enter: executes cell but stays in the same cell
alt+enter: executes cell, insert a new cell below and takes you there
Alternatively, you can use Run
from the menu.
a = 2
b = 6
c = a+b
print(c)
When a cell has been executed, the data are stored and can be used later.
print(a, b)
NB! New cells can be inserted from the menu (Insert
) or by pressing Esc b
. This is useful if you just want to try something without making changes in the existing code.
The code cells are executed in the sequence you choose. There is a number to the left of each cell helping you to keep track of this.
So if you go back and re-execute one of the first cells, data from previous executed cells are still available. This may sometimes cause strange errors. So once in a while it may be useful to restart completely, which can be done from Kernel
in the menu.
Our computations will heavily rely on the modules:
Numpy: Arrays (vectors and matrices), and all the standard mathematical functions.
Matplotlib: Graphs and figures.
SciPy Library: Fundamental library for Scientific computing, including numerical routines for e.g. numerical integration, interpolation, optimization, linear algebra, and statistics.
Occassionally, we might also use
A good entry point to all these libraries is the SciPy web-page.
For our purpose, the first cell will always be something similar to the next one below. All required modules and functions are imported, and Jupyter is made ready for showing plots inside the the document. This cell should always be executed first!
%matplotlib inline
# Import required modules.
import numpy as np
from numpy import pi # Import the number pi
from numpy.linalg import solve # To solve Ax=b
import matplotlib.pyplot as plt
# For plotting.
newparams = {'figure.figsize': (8.0, 4.0), 'axes.grid': True,
'lines.markersize': 8, 'lines.linewidth': 2,
'font.size': 16}
plt.rcParams.update(newparams)
You may change the parameters in newparams
to your own
preferences.
Executing the following cell loads a non-default and more lively css style for the notebook.
Make sure that you download the corresponding css file tma4125.css
from the course web page.
from IPython.core.display import HTML
def css_styling():
styles = open("tma4320.css", "r").read()
return HTML(styles)
# Comment out next line and execute this cell to restore the default notebook style
css_styling()
The implementation of numerical methods is not necessarily very complicated, most of the time we can rely on the following constructions:
Functions
Loops and control statements
Vectors and matrices
Plots and visualization
Numpy is the Python module handling vectors and matrices and all different kind of linear algebra. Vectors and matrices are considered as 1- and 2-dimensional arrays. There is no distinction between row- and column vectors.
The following demonstrates some use of this module. Let $A$ and $b$ be given by
In Python, these are represented by
A = np.array([[ 1.4, 2.2, -1.0], # A: A 2-dimensional array
[ 1.6, -2.7, 1.2],
[-3.2, 1.2, -1.8]])
x = np.array([1.0, -2.0, 3.0]) # x: A 1-dimensional array
print('A = \n', A)
print('\nx = ', x)
Notice that the index starts with 0 in Python, but usually with 1 in mathematics.
Thus $x_3$ in mathematics corresponds to x[2]
in Python, and $a_{21}$ to A[1,0]
.
print(A[1,0])
print(x[2])
If a
and b
are two arrays of the same dimension, standard operations as +
, -
, *
and /
are always element by element operations. The same is the case for a**p
, returning each element of a
taken to the power of p
.
a = np.array([1, 2, 3])
b = np.array([3, 4, 5])
print('a+b = ', a+b)
print('a-b = ', a-b)
print('a/b = ', a/b)
print('a*b = ', a*b)
print('a**b = ', a**b)
print('a**2 = ', a**2)
Similar things happen if you perform this operations between an array
a
and scalar p
, which just means that you perform the
corresponding operation elment-wise. Try it out:
# Insert your code here
The usual matrix-vector product is done by the command A@x
or dot(A,x)
.
y = A@x
print(y)
ip = y@x
print(ip)
And here is a list of quite useful linear algebra functions from the numpy library (start with np.function
)
solve(A,b)
: Solves a linear system $ A \mathbf{x} = \mathbf{b} $. Needs to be imported from numpy.linalg
(which is already done in this notebook)
ones(n)
: returns $[1,1,\dotsc,1]\in\mathbb{R}^n$
zeros(n)
: returns $[0,0,\dotsc,0]\in\mathbb{R}^n$
shape(A)
: The size of an array. Returns $(n, m)$ if $A\in \mathbb{R}^{n\times m}$
eye(n)
: The $ n \times n$ identity matrix $I_n$.
To find
Test them out, and make yourself familiar with them.
# Insert your code here
Contrary to most other programming languages, MATLAB included, there are no end
type expression in Python. Blocks to be executed within a loop or similar have to be indented, and the length of the indentation must be consistent (in these notes we use 4spaces). A block ends when the indentation ends.
In the following, some loop and control structures are demonstrated.
For loop: The aim of this program is to make a nice formatted output of the element of an array.
x = np.array([1.3, 4.6, 2.1, -5.8, 2.3, -3.2]) # The array
n = len(x) # The number of elements in the array
for i in range(n): # Start of the loop
print(f"i = {i:2d}, x = {x[i]:6.2f}") # Formatted output
print('\nAnd this concludes the loop') # End of the loop
Notice that
for i in range(n)
corresponds to $i=0,1,\dotsc,n-1$.
Now write every second term, starting with the second one:
for i in range(1,n,2): # Loops over i=1,3,5
print(f"i = {i:2d}, x = {x[i]:6.2f}") # Formatted output
While loop: Alternatively, you could use a while loop:
i = 0
while (i<n): # Start of the loop
print(f"i = {i:2d}, x = {x[i]:6.2f}") # Formatted output
i += 1 # Increment i by 1
If-statements and break: The loop could be terminated after the second execution:
for i in range(n): # Start of the loop
print(f"i = {i:2d}, x = {x[i]:6.2f}") # Formatted output
if i >= 1: # Break the loop if i >= 1.
break
Let us present an example from physics. Throw a ball straight up in the air. Ignore the air resistance. Let the initial velocity be $v_0$. The position (height) of the ball and its velocity as a function of time $t$ is given by
where $g=9.81\, m/s^2$ is the gravitational accelration at the surface of the earth. The function below takes time, initial velocity and gravitational constant as inputs, and returns the velocity and the position. Notice that the initial velocity and the gravitational constant are given default values.
def ball(t, v0=0, g=9.81):
y = v0*t - 0.5*g*t**2
v = v0 - g*t
return y, v
The function can be used to answer some questions, e.g.:
t = 2.0
y1, v1 = ball(t)
print('v0=',0.0,', t=', t, ', y=', y1, ', v=', v1)
v_start = 1.0
y2, v2 = ball(t, v0=v_start)
print('v0=', v_start,', t=', t, ', y=', y2, ', v=', v2)
y3, v3 = ball(t, v0=v_start, g=1.625)
print('v0=', 1.0,', t=', t, ', y=', y3, ', v=', v3, ' on the moon.')
Whenever possible, write your functions such that if the input is an array, so is the output. This is the case here, and it makes it easy to generate a graph of the function, e.g., the position as a function of time when the initial velocity is 10 m/s.
t = np.linspace(0,2,101)
v_start = 10
y, v = ball(t, v0=v_start)
plt.plot(t,y)
plt.xlabel('t')
plt.ylabel('y')
plt.title('Position of the ball when v0 = {} m/s'.format(v_start));
This note provides only a very short introduction to some useful functions and constructions in Python. More will be used throughout the course. We hope that most of it will be understandable or explained along the course, but if there are functions or structures you do not understand, please look them up yourself, for instance by
?function
.?np.linspace
If you want to see what functions are available in an already imported library, e.g. Numpy and Matplotlib, write
np.<TAB>
plt.<TAB>
and you will get more information than you want. It is probably much more efficient to read the documentation:
Time spent on getting acquainted with the documentation of these libraries is very well invested time, not only for this course.
For a nice introduction for scientific computing in Pyhton, we can recommend