- \n",
"
- Overview and Learning Objectives \n", "
- NumPy Module \n", "
- NumPy Arrays \n", "
- Array Reshaping \n", "
- Array Operations \n", "
- NumPy Type Conversion \n", "
- NumPy Constants \n", "
- Exercise 1: NumPy Array Manipulations \n", "
- Recap: Matrices and Basic Operations \n", "
- Exercise 2: Matrix Operations \n", "
- Exercise 3: Mathematical Functions \n", "

\n",
"## Overview and Learning Objectives

\n",
"\n",
"Python has several useful built-in packages as well as additional external packages. One such package is NumPy, which adds support for multi-dimensional arrays and matrices, along with a number of mathematical functions to operate on these structures. This unit covers array objects as the most fundamental NumPy data structure along with important NumPy array operations. Furthermore, we discuss numerical types, methods for type conversion, and constants (including the constants `nan` and `inf`) offered by NumPy. The three exercises included at the end of this unit cover key aspects needed throughout the PCP notebooks. Therefore, we encourage students to work through these exercises carefully. In particular, we recapitulate the mathematical concept of matrices and matrix multiplication, before we cover these aspects in Exercise 2 from an implementation perspective. We believe that understanding both—the mathematical formulation of abstract concepts and their realization in a programming context—is key in engineering education. This philosophy also forms the basis of the PCP notebooks to follow. \n",
" \n",
"

"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"## NumPy Module\n",
"\n",
"As said above, NumPy is a Python library used for working with arrays. In principle, one can also use the Python concept of lists to model an array. However, processing such lists is usually slow. Being mostly written in C or C++ and storing arrays at continuous places in memory (unlike lists), NumPy can process and compute with arrays in a very efficient way. This is the reason why NumPy is so powerful. In order to use NumPy, one needs to install the `numpy` package. This is what we already did when we created the [PCP Environment](PCP_01_getstarted.html), which contains a `numpy`-version. Furthermore, we need to use the Python `import` statement to get access to the NumPy functions. It is convenient to assign a short name to a frequently-used package (for example `np` in the case of `numpy`). This short name is used as a prefix when calling a function from the package. In the following code cell, we import `numpy` as `np`. To find out what the `numpy`-module contains, we use the Python command `dir(np)` to return a list all properties and methods of the specified module."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"list_numpy_names = dir(np)\n",
"print(list_numpy_names)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"## NumPy Arrays\n",
"\n",
"The fundamental NumPy data structure is an **array** object, which represents a multidimensional, homogeneous array of fixed-size items. Each array has a shape, a type, and a dimension. One way to create a NumPy array is to use the `array`-function. In the following code cells, we provide various examples and discuss basic properties of NumPy arrays. "
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"\n",
"x = np.array([1, 2, 3, 3])\n",
"print('x = ', x)\n",
"print('Shape:', x.shape)\n",
"print('Type:', x.dtype)\n",
"print('Dimension:', x.ndim)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Note that, in this example, `x.shape` produces a one-element tuple, which is encoded by `(4,)` for disambiguation. (The object `(4)` would be an integer of type `int` rather than a tuple.) Two-dimensional arrays (also called **matrices**) are created like follows:"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"x = np.array([[1, 2, 33], [44, 5, 6]])\n",
"print('x = ', x, sep='\\n')\n",
"print('Shape:', x.shape)\n",
"print('Type:', x.dtype)\n",
"print('Dimension:', x.ndim)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"There are a couple of NumPy functions for creating arrays:"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"print('Array of given shape and type, filled with zeros: ', np.zeros(2))\n",
"print('Array of given shape and type, filled with integer zeros: ', np.zeros(2, dtype='int'))\n",
"print('Array of given shape and type, filled with ones: ', np.ones((2, 3)), sep='\\n')\n",
"print('Evenly spaced values within a given interval: ', np.arange(2, 8, 2))\n",
"print('Random values in a given shape: ', np.random.rand(2, 3), sep='\\n')\n",
"print('Identity matrix: ', np.eye(3), sep='\\n')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"## Array Reshaping \n",
"\n",
"Keeping the total number of entries, there are various ways for reshaping an array. Here are some examples:"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"x = np.arange(2 * 3 * 4)\n",
"print('x =', x)\n",
"print('Shape:', x.shape)\n",
"\n",
"y = x.reshape((3, 8))\n",
"print('y = ', y, sep='\\n')\n",
"print('Shape:', y.shape)\n",
"\n",
"z = np.reshape(x, (3, 2, 4))\n",
"print('z = ', z, sep='\\n')\n",
"print('Shape:', z.shape)\n",
"\n",
"print('Element z[0, 1, 2] = ', z[0, 1, 2])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"NumPy allows for giving one of the new shape parameters as `-1`. In this case, NumPy automatically figures out the unknown dimension. Note the difference between the shape `(6,)` and the shape `(6,1)` in the following example."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"x = np.arange(6)\n",
"print(f'Shape: {x.shape}; dim: {x.ndim}')\n",
"x = x.reshape(-1, 2)\n",
"print(f'Shape: {x.shape}; dim: {x.ndim}')\n",
"x = x.reshape(-1, 1)\n",
"print(f'Shape: {x.shape}; dim: {x.ndim}')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"## Array Operations \n",
"\n",
"There are many ways to compute with NumPy arrays. Many operations look similar to the ones when computing with numbers. Applied to arrays, these operations are conducted in an element-wise fashion:"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"x = np.arange(5)\n",
"print('x = ', x)\n",
"print('x + 1 =', x + 1)\n",
"print('x * 2 =', x * 2)\n",
"print('x * x =', x * x)\n",
"print('x ** 3 =', x ** 3)\n",
"print('x / 4 =', x / 4)\n",
"print('x // 4 =', x // 4)\n",
"print('x > 2 =', x > 2)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Here are some examples for computing with matrices of the same shape. It is important to understand the difference between element-wise multiplication (using the operator `*`) and usual matrix multiplication (using the function `np.dot`):"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"a = np.arange(0, 4).reshape((2, 2))\n",
"b = 2 * np.ones((2, 2))\n",
"print('a = ', a, sep='\\n')\n",
"print('b = ', b, sep='\\n')\n",
"print('a + b = ', a + b, sep='\\n')\n",
"print('a * b (element-wise multiplication) = ', a * b, sep='\\n')\n",
"print('np.dot(a, b) (matrix multiplication) = ', np.dot(a, b), sep='\\n')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Note that arrays and lists may behave in a completely different way. For example, using the operator `+` leads to the following results:"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"a = np.arange(4)\n",
"b = np.arange(4)\n",
"print(a + b, type(a + b))\n",
"\n",
"a = list(a)\n",
"b = list(b)\n",
"print(a + b, type(a + b))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The sum of an array's elements can be computed along the different dimensions specified by the parameter `axis`. This is illustrated by the following example:"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"x = np.arange(6).reshape((2, 3))\n",
"print('x = ', x, sep='\\n')\n",
"print('Total sum:', x.sum())\n",
"print('Column sum: ', x.sum(axis=0))\n",
"print('Row sum:', x.sum(axis=1))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"There are many ways for accessing and manipulating arrays. Here are some examples:"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"x = np.arange(6).reshape((2, 3))\n",
"print('x = ', x, sep='\\n')\n",
"print('Element in second column of second row:', x[1, 1])\n",
"print('Boolean encoding of element positions with values larger than 1: ', x > 1, sep='\\n')\n",
"print('All elements larger than 1:', x[x > 1])\n",
"print('Second row:', x[1, :])\n",
"print('Second column:', x[:, 1])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"## NumPy Type Conversion\n",
"\n",
"In the [PCP Notebook on Python Basics](PCP_02_python.html), we have already discussed the standard numeric types `int`, `float`, and `complex` offered by Python (and the function `type()` to identify the type of a variable). The NumPy package offers many more numerical [types and methods for type conversion](https://numpy.org/doc/stable/user/basics.types.html). We have already seen how to obtain the type of a numpy array using `dtype`. One can create or convert a variable to a specific type using `astype`. Some examples can be found in the next code cell."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"a = np.array([-1,2,3])\n",
"print(f'a = {a}; dtype: {a.dtype}')\n",
"b = a.astype(np.float64)\n",
"print(f'b = {b}; dtype: {b.dtype}')\n",
"c = a.astype(np.int64)\n",
"print(f'c = {c}; dtype: {c.dtype}')\n",
"d = a.astype(np.uint8)\n",
"print(f'd = {d}; dtype: {d.dtype}')\n",
"e = np.array([1,2,3]).astype(np.complex64)\n",
"print(f'e = {e}; dtype: {e.dtype}')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Specifying the exact type is often important when using packages such as [`numba` for optimizing machine code at runtime](http://numba.pydata.org/). In the following example, we give an example where a wrong initialization leads to an error (or warning with some `nan` results) when computing the square root of a negative number."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"print('=== Initialization with \\'int32\\' leading to an error ===', flush=True)\n",
"x = np.arange(-2, 2)\n",
"print(x, x.dtype)\n",
"x = np.sqrt(x)\n",
"print(x) \n",
"\n",
"print('=== Initialization with \\'complex\\' ===', flush=True)\n",
"x = np.arange(-3, 3, dtype='complex')\n",
"print(x, x.dtype)\n",
"x = np.sqrt(x)\n",
"print(x)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"## NumPy Constants\n",
"\n",
"NumPy offers several constants that are convenient to compute with. In the following, we give some examples:"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"print(f'Archimedes constant Pi: {np.pi}; type: {type(np.pi)}')\n",
"print(f'Eulerâ€™s constant, base of natural logarithms: {np.e}; type: {type(np.e)}')\n",
"print(f'Floating point representation of (positive) infinity: {np.inf}; type: {type(np.inf)}')\n",
"print(f'Floating point representation of (negative) infinity: {np.NINF}; type: {type(np.NINF)}')\n",
"print(f'floating point representation of Not a Number (NaN): {np.nan}; type: {type(np.nan)}')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"In particular, the constants `nan` and `inf` can be convenient to avoid case distinctions. However, computing with such constants can also be a bit tricky as the following examples show."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"a = 10\n",
"b = np.inf\n",
"c = -np.inf\n",
"d = np.nan\n",
"\n",
"print(f'a = {a}, b = {b}, c = {c}, d = {d}')\n",
"print('a + b =', a + b)\n",
"print('a * b =', a * b)\n",
"print('a + c =', a + c)\n",
"print('a - c =', a - c)\n",
"print('a + d =', a + d)\n",
"print('b + c =', b + c)\n",
"print('b * c =', b * c)\n",
"print('b / c =', b / c)\n",
"print('b + d =', b + d)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
" NumPy offers functions such as `np.where`

and `np.isinf`

to check for special constants. In the following, the class [`np.errstate`](https://numpy.org/doc/stable/reference/generated/numpy.errstate.html) is used to suppress warning that are usually output when dividing by zero. "
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"print('Test element-wise for positive or negative infinity:', np.isinf([np.inf, np.NINF, np.nan]))\n",
"\n",
"a = np.arange(-2, 3) \n",
"print('a = ', a)\n",
"with np.errstate(divide='ignore', invalid='ignore'):\n",
" b = a / 0\n",
"print('b = a / 0 =', b)\n",
"\n",
"ind_inf = np.isinf(b)\n",
"ind_inf_pos = np.where(ind_inf)\n",
"print('Indices with inf values:', ind_inf, ind_inf_pos)\n",
"\n",
"ind_nan = np.isnan(b)\n",
"ind_nan_pos = np.where(ind_nan)\n",
"print('Indices with nan values:', ind_nan, ind_nan_pos)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Exercises and Results"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import libpcp.numpy\n",
"show_result = True"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"**Exercise 1: NumPy Array Manipulations**

\n", "

"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#\n", "

- \n",
"
- Create a NumPy array
`a`

with ascending natural numbers in the interval $[10, 20]=\\{10,11,\\ldots,20\\}$ (using`np.arange`

). \n",
" - Set all entries of
`a`

to zero where`a`

$\\leq13$ and`a`

$>16$. \n",
" - Extend the resulting array
`a`

with a NumPy array containing the numbers of the interval $[4,6]$ and store the result in a variable`b`

(using`np.append`

). \n",
" - Remove duplicate values of
`b`

(using`np.unique`

) and store the result in a variable`c`

. Note that the result is automatically sorted in ascending order. \n",
" - Sort
`c`

in descending order and store the result in a variable`d`

. Explore and discuss various options to do this including the slicing method`c[::-1]`

, the function`reversed()`

, as well as the NumPy functions`np.sort`

and`np.flip`

. \n",
"

\n",
"**Recap: Matrices and Basic Operations**

\n", "Let $N,M\\in\\mathbb{N}$ be two positive integers. An $(N\\times M)$-matrix $\\mathbf{A}$ is a rectangular array of entries arranged in**rows** and **columns**. For example, if entries are real numbers, one also writes $\\mathbf{A}\\in\\mathbb{R}^{N\\times M}$. Let $a_{nm}\\in\\mathbb{R}$ denote the entry of $\\mathbf{A}$ being in the $n^{\\mathrm{th}}$ row and $m^{\\mathrm{th}}$ column for $n\\in[1:N]=\\{1,2,...,N\\}$ and $m\\in[1:M]$. Then, one also often writes $\\mathbf{A}=[a_{nm}]$. A **vector** is a matrix where either $N=1$ (then called **row vector**) or $M=1$ (then called **column vector**).\n",
" \n",
"When computing with vectors and matrices, one has to pay attention that the dimensions $N$ and $M$ of the operands match properly. For example, for a **matrix summation** or **matrix subtraction**, the operands must be of same dimensions or one of them needs to be a scalar (in this case the operations are applied per entry). The multiplication of matrices refers to a **matrix multiplication**. For an $(N\\times M)$-matrix $\\mathbf{A} = [a_{nm}]$ and a $(M\\times P)$-matrix $\\mathbf{B} = [b_{mp}]$, the product matrix $\\mathbf{C} = \\mathbf{A}\\mathbf{B}$ with entries $\\mathbf{C}=[c_{np}]$ is defined by $c_{np} = \\sum_{m=1}^M a_{nm}b_{mp}$ for $n\\in[1:N]$ and $p\\in[1:P]$. In other words, the entry $c_{np}$ is the **inner product** (sometimes also called **scalar product**) of $n^{\\mathrm{th}}$ row of $\\mathbf{A}$ and $p^{\\mathrm{th}}$ column of $\\mathbf{B}$. This calculation rule is illustrated by the following figure:\n",
"\n",
"\n",
"

\n", "

"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n", "Let $N,M\\in\\mathbb{N}$ be two positive integers. An $(N\\times M)$-matrix $\\mathbf{A}$ is a rectangular array of entries arranged in

\n", "

\n",
"**Exercise 2: Matrix Operations**

\n", "

"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#\n", "

- \n",
"
- Construct a matrix $\\mathbf{B} = \\begin{bmatrix}2 & 2 & 2 & 2\\\\2 & 2 & 2 & 2\\\\0 & 0 & 0 & 0\\\\\\end{bmatrix}$ just using the NumPy functions
`np.zeros`

,`np.ones`

, and`np.vstack`

. Check the matrix dimensions using`np.shape`

. \n",
" - Find the row and column index of the maximum entry of the matrix $\\mathbf{D} = \\begin{bmatrix}2 & 0 & 2\\\\-1 & 5 & 10\\\\-3 & 0 & 9\\\\\\end{bmatrix}$ using the functions
`np.max`

,`np.argmax`

and`np.unravel_index`

. Why is it not sufficient to use`np.argmax`

? \n",
" - Given a row vector $\\mathbf{v} = [3\\;2\\;1]$ and a column vector $\\mathbf{w} = [6\\;5\\;4]^T$, compute $\\mathbf{v}\\mathbf{w}$ and $\\mathbf{w}\\mathbf{v}$ using
`np.dot`

and`np.outer`

. What is the difference between`np.multiply`

,`np.dot`

, and`np.outer`

? \n",
" - Given $\\mathbf{A} = \\begin{bmatrix}1 & 2\\\\3 & 5\\end{bmatrix}$, $\\mathbf{v} = \\begin{bmatrix}1\\\\4\\end{bmatrix}$, compute $\\mathbf{A}^{-1}$ and $\\mathbf{A}^{-1}\\mathbf{v}$ (using
`np.linalg.inv`

). \n",
"

\n",
"**Exercise 3: Mathematical NumPy Functions**

\n", "The NumPy package offers many different mathematical functions. Explore these functions by trying them out on small examples. In particular, you should gain a good understanding of the following concepts.\n", " \n", "

"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#\n", "The NumPy package offers many different mathematical functions. Explore these functions by trying them out on small examples. In particular, you should gain a good understanding of the following concepts.\n", " \n", "

- \n",
"
- Generate a NumPy array
`v_deg`

containing the numbers $[0, 30, 45, 60, 90, 180]$. Apply the function`np.deg2rad`

to convert this array (given in degree) to an array`v_rad`

(given in radiants). Then apply the functions`np.cos`

and`np.sin`

and discuss the results. \n",
" - Using the the same array as before, apply the exponential function
`np.exp(1j * v_rad)`

. What is meant by`1j`

? What is the relation to`np.cos`

and`np.sin`

? Use the functions`np.real`

,`np.imag`

, and`np.isclose`

to make this relation explicit.\n", " - Try out different functions for rounding using the numbers $[-3.1416, -1.5708, 0, 1.5708, 3.1416]$. What is the difference between the functions
`np.round`

,`np.floor`

,`np.ceil`

, and`np.trunc`

? \n",
"