0% found this document useful (0 votes)
8 views165 pages

Python For Engineers

The document titled 'Python for Engineers' by Lorenz Kies provides a comprehensive guide to programming in Python, covering fundamental concepts such as variables, functions, loops, and data structures. It also includes sections on advanced topics like NumPy and Matplotlib for numerical and graphical data analysis. The content is structured to facilitate learning for engineers and includes practical examples and exercises.

Uploaded by

iampaulisaac
Copyright
© All Rights Reserved
We take content rights seriously. If you suspect this is your content, claim it here.
Available Formats
Download as PDF, TXT or read online on Scribd
0% found this document useful (0 votes)
8 views165 pages

Python For Engineers

The document titled 'Python for Engineers' by Lorenz Kies provides a comprehensive guide to programming in Python, covering fundamental concepts such as variables, functions, loops, and data structures. It also includes sections on advanced topics like NumPy and Matplotlib for numerical and graphical data analysis. The content is structured to facilitate learning for engineers and includes practical examples and exercises.

Uploaded by

iampaulisaac
Copyright
© All Rights Reserved
We take content rights seriously. If you suspect this is your content, claim it here.
Available Formats
Download as PDF, TXT or read online on Scribd

Python for Engineers

Lorenz Kies

Apr 13, 2026


CONTENTS

1 Getting Started 9
1.1 My first Program: “Hello World” . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
1.2 A Calculator on Steroids . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
1.3 Variables & Expressions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
1.4 Types & Values . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12
1.5 Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
1.6 Conditionals with if and else . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
1.7 Operators . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
1.7.1 Arithmetical Operators . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
1.7.2 Relational Operators . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
1.7.3 Logical Operators . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
1.7.4 Operator Precedence . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
1.8 Examples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
1.8.1 Beam Deflection . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
1.8.2 Maximum . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
1.8.3 Price Calculation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21
1.9 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21

2 Moving On 22
2.1 Loops . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22
2.1.1 for Loops . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22
2.1.2 The range Function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23
2.1.3 Iterating Over a Collection . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
2.1.4 while Loops . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
2.1.5 Executing Code Step by Step with a Debugger . . . . . . . . . . . . . . . . 26
2.1.6 break and continue . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
2.2 More About Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
2.2.1 Pure Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
2.2.2 Docstrings . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29
2.2.3 Scopes and Variable Visibility . . . . . . . . . . . . . . . . . . . . . . . . . 30
2.3 Examples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
2.3.1 Factorial . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32
2.3.2 Herons Algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
2.4 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34

3 Built-in Data Structures 35


3.1 Lists . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35
3.1.1 Operations on Lists . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37
3.1.2 Multi-Dimensional Lists . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39

i
3.1.3 List Examples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40
3.2 Tuples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
3.3 Dictionaries . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43
3.4 Sets . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46
3.5 More Examples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47
3.5.1 Why Sets over Lists? . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47
3.6 Dictionaries for Caching . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
3.7 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49

4 NumPy & Matplotlib 50


4.1 Installing Packages . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50
4.2 NumPy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
4.2.1 Note on Dot Notation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
4.2.2 NumPy Arrays . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52
4.2.3 Vectorized Computations . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53
4.2.4 Array Creation Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55
4.2.5 Linear Algebra . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56
4.2.6 Examples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
4.3 Matplotlib . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
4.3.1 Basic Plotting . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
4.3.2 Subplots . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64
4.3.3 Examples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65
4.4 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67

5 Strings, Files & File Systems 68


5.1 Task Description . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68
5.2 Working with the File System . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69
5.3 Reading and Writing Files . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73
5.4 Working with Strings . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74
5.5 String Interpolation with f-Strings . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76
5.6 Putting it all together . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 77
5.7 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79

6 Object-Oriented Programming 80
6.1 Defining a Class . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80
6.2 Special Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81
6.3 Complex Numbers in Python . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82
6.4 Encapsulation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83
6.5 Inheritance . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 84
6.6 Dataclasses & Namedtuples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87
6.7 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 88
6.8 Appendix - Decorators . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 88

7 Bits and Pieces 90


7.1 Unpacking Assignments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90
7.2 Better Iteration with zip() and enumerate() . . . . . . . . . . . . . . . . . . . . 91
7.3 List Comprehensions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92
7.4 Reference Semantics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 93
7.5 Lambdas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96
7.6 Function Parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97
7.7 Code Style . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99
7.7.1 Naming Conventions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99

ii
7.7.2 Automatic Formatting & Linting . . . . . . . . . . . . . . . . . . . . . . . . . 100
7.7.3 Use More Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100
7.7.4 Use Early Returns . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100
7.7.5 Use Type Hints . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101
7.7.6 Prioritize Readability . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101
7.7.7 Don’t Be Afraid to Start Over . . . . . . . . . . . . . . . . . . . . . . . . . . 101
7.8 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101

8 Pandas 103
8.1 Series for Ordered Data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103
8.2 DataFrames as Collections of Series . . . . . . . . . . . . . . . . . . . . . . . . 104
8.3 Convenience Functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 107
8.4 Indexing & Filtering . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 107
8.5 Saving & Loading DataFrames . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 108
8.6 Grouping & Aggregating . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 109
8.7 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 113

9 Advanced Plotting 114


9.1 Plots for Visualizing Distributions of Data . . . . . . . . . . . . . . . . . . . . . . . . 114
9.2 Interactive Plots with Jupyter Widgets . . . . . . . . . . . . . . . . . . . . . . . . . 116
9.3 Matplotlib Interactive Backend . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117
9.4 3D Plotting . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119
9.5 Animations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 121
9.6 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123

10 Differential Equations & Systems 124


10.1 First Order Differential Equations . . . . . . . . . . . . . . . . . . . . . . . . . . . . 125
10.2 Higher Order Differential Equations . . . . . . . . . . . . . . . . . . . . . . . . . . . 126
10.3 Differential Equations with Inputs . . . . . . . . . . . . . . . . . . . . . . . . . . . . 127
10.4 Linear Time-Invariant Systems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 128
10.4.1 State Space Representation . . . . . . . . . . . . . . . . . . . . . . . . . . 129
10.4.2 Transfer Function Representation . . . . . . . . . . . . . . . . . . . . . . . 132
10.4.3 Other Operations on LTI Systems . . . . . . . . . . . . . . . . . . . . . . . 132
10.5 Phase Portraits . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 134
10.6 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 136

11 Root Finding & Optimization 137


11.1 Solving Equations (Root Finding) . . . . . . . . . . . . . . . . . . . . . . . . . . . . 137
11.1.1 Scalar Problems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 137
11.1.2 Limitations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 139
11.1.3 Multiple Unknowns and Equations . . . . . . . . . . . . . . . . . . . . . . . 140
11.2 Optimization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 141
11.3 Curve Fitting . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 143
11.4 Optimal Control . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 145
11.4.1 Offline Optimization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 145
11.4.2 Online Optimization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 148
11.5 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 149

12 Python Projects 150


12.1 Importing Code from Other Files . . . . . . . . . . . . . . . . . . . . . . . . . . . . 150
12.1.1 Files in the Same Directory . . . . . . . . . . . . . . . . . . . . . . . . . . . 150
12.1.2 Files in Different Directories . . . . . . . . . . . . . . . . . . . . . . . . . . . 152

iii
12.2 Making Scripts Importable . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 155
12.3 Minimal Project Setup . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 156
12.4 Managing Dependencies . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 157
12.5 Testing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 157
12.6 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 160

iv
Python for Engineers

Lorenz Kies digitale-werkzeuge@[Link]


Links:
• Book PDF
• Book Website
• Solved Lectures
• Blank Lectures
• Exercises
• JupyterLite
• eduVote

Introduction

This course is meant as an introduction to Python and programming in general for engineers. We
live in an age where practically unimaginable computing resources are almost universally avail-
able. Being able to make use of these resources has become a requirement to work productively
in many science and engineering related fields. In addition to that, programming might also help
you in other avenues of your life. It is therefore important, that you as engineers have a basic
understanding of programming and some commonly available tools.
Python as a language and ecosystem comes with a lot of libraries, tools and other resources for
almost everything you might want to do. In addition to that it is a fairly beginner-friendly language,
which make this is a great candidate for this course.
In this course you will learn:
• The fundamentals of programming in Python (variables and expressions, control structures,
functions, common data structures)
• The basics Python’s rich scientific tool set (NumPy, SciPy, Matplotlib, Pandas, Python Con-
trol)
• How to find and read documentation
• How to use the tools required for developing and running Python code
• And how to use those skills to solve everyday engineering problems
Additionally, there will be a lot of pointers for further reading and learning to lead you to more
advanced topics.

Course Structure

The course is structured into:


• 12 live coding lectures & tutorials
• 4 exercise sheets and 1 mini project (solutions presented in tutorials)
• 2 graded programming tests (probably on 04.06.2026 and 16.07.2026 14:00)

CONTENTS 1
Python for Engineers

• And multiple opportunities to work with real world experiments


The course as a whole is on pass/fail basis you will not get a number grade. To pass you need to
achieve at least 50% of the respective points in across the two graded tests.
The lectures are held as live coding sessions in Jupyter notebooks where you the students should
type along. Questions and other forms of interaction are welcome. Don’t hesitate to ask curious
questions beyond the immediate scope of the lecture, sidetracks are welcome and should make
for a more lively learning experience.
The exercise sheets are also (mostly) Jupyter notebooks and are meant to be solved before the
respective tutorial. It is recommended, that you work on these in pairs. You can use tools like
VS Code Live Share to work together remotely if needed. The mini project should give more
insight into Python outside of Jupyter notebooks and should also be done in pairs. You should
consider these exercises as mandatory, unlike some other courses, doing nothing until
a few days before the test is absolutely not a feasible strategy. Learning to program
requires consistent practice.
There will be additional bonus exercises in the form of “Hardware Challenges”, if you submit
working code for these problems, we will invite you to the lab where you can run your code on
the actual experiment.

Schedule (Preliminary)

Date Lecture Tutorial


14.04.2026 Getting Started -
21.04.2026 Moving On Exercise 0.1 to 0.3
28.04.2026 Built-in Data Structures Exercise 0.4 to 0.7
05.05.2026 NumPy & Matplotlib Exercise 1.1.1 to 1.1.9
12.05.2026 Strings, Files & File Systems Exercise 1.2.1 to 1.4.1
19.05.2026 Object-Oriented Programming Exercise 2.1.1 to 2.2.2
26.05.2026 - -
02.06.2026 Bits and Pieces Test Preparation
09.06.2026 Pandas Exercise 2.4 to 2.5.2
16.06.2026 Advanced Plotting Exercise 2.6.1 & 2.7.1 & 3.1
23.06.2026 Differential Equations & Systems Exercise 3.2.1 & 3.2.2
30.06.2026 Root Finding & Optimization Exercise 3.3.1 to 3.3.5
07.07.2026 Python Projects Exercise 3.4.1 to 3.4.4
14.07.2026 Buffer Exercise 4.1 to 4.7

Why Is All of This in English?

All of the course material is in English, this is not to inconvenience you, but it should ideally be
to help you. Pretty much all programming languages as well as related tools and most online
resources are in English. This makes programming as a whole mostly an exclusively English
thing. To help you get used to programming in English i.e. English variable names, interfaces and

CONTENTS 2
Python for Engineers

documentation the course material is in English. To not pose too much of an additional challenge
to people who are not as comfortable in English, the lectures will be held in German.

AI Tools - A Word of Caution

Large language model based AI tools have seen a rapid increase in capabilities and popularity in
recent years. This is very much a double-edged sword in the context of learning to program. On
the one hand they allow you to do more things faster with less knowledge, but on the other hand
they allow you to do more things faster with less knowledge. In a few more words:
AI tools are now capable of solving many small common problems decently, so they can be an at-
tractive option if you are not interested in learning and just need quick results. They are especially
good at solving most of the classical introduction to programming exercises. This might easily
you lead to believe, that there is no use in learning “how to generate a list of Fibonacci numbers”
since you could always just have an AI do it, just like you don’t do arithmetic with pen and paper
since you have a calculator. But solving these problems is not about having generating a solution
to these problems, but about the skills you learn along the way. Skills that you will need to solve
more complex problems even with the help of AI tools.
Here is an analogy to illustrate this idea: In your very first swimming lessons, you will likely be in
a pool where you can easily stand up (i.e. use AI tools). So you might ask why you would need to
learn to swim when you can just stand up. But you are not learning to swim so you can keep you
head above water in the beginners pool (i.e. solve practice problems) you are learning to swim so
you won’t drown when you can’t stand (i.e. solve real world problems).
Not only are these tools very good at solving practice problems, they can be even worse: depend-
ing on how tightly they are integrated into your workflow, they might not only deprive you of the
valuable learning experience, but fool you into believing you learned something that you didn’t. So
while AI tools, if used right, can be a helpful tool in learning to program, they can quickly become
a crutch that hinders learning. That is why, for learning to program, I recommend that you
do not use any AI tools.

Setting up your Environment

To be able to follow the live lectures and solve the exercise sheets you will need to set up your de-
velopment environment. It is recommended that you install Python using uv and use Visual Studio
Code as your code editor. The installation of those tools is described in the following sections.

Shells & Package Managers

For the installation of programming environments, but also for programming itself, it is very useful
to have at least some understanding of shells / command line interfaces and package managers.
Shells are text-based interfaces to the operating system, which allow you to interact with the sys-
tem by typing commands into a terminal and getting text-based output back.
While there are also other options, the shells you will typically be using are:
• PowerShell on Windows

CONTENTS 3
Python for Engineers

• Bash on Linux (the application is usually called “Terminal”)


• Zsh on MacOS (the application is also called “Terminal”)
You can launch them using the search function of your operating system.
In these terminals you can then use package mangers to install applications and other software.
Package managers are kind of like the app stores on your phones only that they are text-based
and much more comprehensive, but they mostly fulfill mostly the same purpose.
Again, there are many alternatives but the common ones are:
• WinGet for Windows
• APT for Linux
• Homebrew for MacOS (brew for short)
Package managers are not limited to installation of applications for operating systems, they are
also used for programming languages. But unlike package managers for operating systems, these
usually don’t install applications, they mostly install libraries which are basically common pieces
of code that you wouldn’t want to rewrite yourself every time. For Python there are PIP, Conda
and a few more.
In the following you should install the required applications using the package managers men-
tioned above.

Setting up Python with uv

There are many ways to install Python on your system. As mentioned above, we will be using uv
which includes an alternative implementation of PIP in addition to some more advanced features
for project management.

Windows

Install uv using WinGet in PowerShell:

winget install [Link]

You may have to agree to the terms of service the first time you use WinGet.

Linux & MacOS

Download and execute the installation script in the terminal:

curl -LsSf [Link] | sh

or using wget if you don’t have curl:

wget -qO- [Link] | sh

CONTENTS 4
Python for Engineers

Create the Python Environment

First, test if your installation of uv was successful by entering the following command in your
terminal:

uv

You should see an output similar to this:

An extremely fast Python package manager.


Usage: uv [OPTIONS] <COMMAND>
...

Unlike many other methods of installing Python, uv does not install Python globally on your system.
Instead, uv is used to create isolated Python installations called virtual environments for each
project. For these purposes, we will consider this course to be the project. Start by creating
a directory where you want to store the material for this course. Then open your shell in said
directory. You can usually do this by right-clicking in your file explorer while holding the shift key
and selecting “Open PowerShell window here” or “Open in Terminal” or something similar.
Start by verifying that you are in the correct directory with:

pwd

This should display the path to the directory you just created.
If the location is correct, initialize the project with uv:

uv init . --bare --name python-for-engineers

Then add the required packages to the project:

uv add numpy scipy matplotlib pandas control ipykernel ipywidgets ipympl␣


↪pyarrow otter-grader pytest ruff vscdark

This will add them to your project specification in stored in a file called [Link] and
install them in a virtual environment for this project.

Executing Python Code

In the first few lectures, we will only be using Jupyter notebooks, so you can skip this section for
now.
There are two primary ways to execute Python code in this environment.
The first option is to execute Python through uv:

uv run python <your_script.py>

Alternatively, you can activate the virtual environment created by uv and then execute Python
directly:
Using Windows PowerShell:

CONTENTS 5
Python for Engineers

.\venv\Scripts\Activate.ps1

You may first need to allow the execution of scripts by opening a new PowerShell window as
administrator and running:

Set-ExecutionPolicy RemoteSigned

On Linux and MacOS you can activate the environment with:

source venv/bin/activate

Setting up Visual Studio Code

To edit your code, you will need a code editor. For this we will use Visual Studio Code. You can
also use the graphical installer from the Visual Studio Code website if you are having trouble with
the instructions.

Windows

Install Visual Studio Code using WinGet in PowerShell:

winget install vscode

Linux

Install with snap in the terminal:

sudo snap install --classic code

MacOS

Install with brew in the terminal:

brew install --cask visual-studio-code

Extensions

So set up VS-Code for Python development you will need to install some extensions. You can
open the extensions tab by pressing Ctrl+Shift+X or by clicking the respective icon in the
sidebar.

CONTENTS 6
Python for Engineers

For the course you will need:


• [Link] Python
• [Link] Jupyter
Install them by clicking the blue “Install” button next to the extension.
It is also recommended to install the following extensions:
• [Link] Ruff
• [Link] Error Lens
• [Link] autoDocstring
• [Link]-spell-checker Code Spell Checker
• [Link]-all-in-one Markdown All in One
• [Link]-csv Rainbow CSV

Using Visual Studio Code

To use all of VS Code’s features, it needs to know what files you are working on and what Python
environment you are using. To do this, you need to open the folder you created for this course as
a workspace in VS Code. You can do this by clicking File -> Open Folder... and selecting
the folder you created for this course.
Next you need to make sure, that VS Code is using the correct Python environment. When you
open a Python file, you should see a notification in the bottom right corner asking you to select
a Python environment. Click on it and select the environment that is located in the venv folder
in your project directory. If you don’t see the notification, you can also select the environment by
clicking on the Python version displayed in the bottom left corner. Alternatively, you can also open
the dialog for selecting the Python environment with the command palette (Ctrl+Shift+P) and
searching for “Python: Select Interpreter”.
Similarly, when you open a Jupyter notebook, you should see a notification asking you to select a
kernel. Click on it and select the kernel that is located in the venv folder in your project directory.

CONTENTS 7
Python for Engineers

If you don’t see the notification, you can also select the kernel by clicking on the kernel name
displayed in the top right corner of the notebook editor.

Shortcuts

Code editors like VS-Code come with a lot of shortcuts to make editing code faster and less te-
dious. Here are some of the most useful ones, you can find a more comprehensive list here.

Key Combination Action


Ctrl+←/→ Move cursor left/right by one word
Pos1/Ende Move cursor to beginning/end of line
Ctrl+# Comment / uncomment line or selection
Alt+Shift+↑/↓ Duplicate line or selection up/down
Alt+↑/↓ Move line or selection up/down
Ctrl+D Add cursor to next occurrence of selected text
F2 Rename symbol (variable, function, class, etc.)
Ctrl+Enter Execute current cell in Jupyter notebook
Shift+Enter Execute current cell and move to next cell in Jupyter notebook
Ctrl+Shift+P Open command palette

CONTENTS 8
CHAPTER

ONE

GETTING STARTED

Open in JupyterLite: Solved Blank


It is time to get started programming in Python. In this first lecture we will get to know Python and
how to use it as a more advanced calculator and seeing how some of Python’s features relate to
concepts from mathematics that are already familiar with.

1.1 My first Program: “Hello World”

A common first program that to write when learning a new language is the “Hello World” program.
It simply prints out the text "Hello World" and in Python it is only one line. And print, in this
case, does not mean to literally print something out (not anymore), but to display the text on the
standard output, meaning the screen. We can execute a cell by pressing Ctrl + Enter or by
clicking the play button in top left of the cell.

print("Hello World")

Hello World

1.2 A Calculator on Steroids

In some sense, you could think of Python as a kind of like a calculator on steroids. It can do
most things that you are used to from normal calculators. Like basic arithmetic and other com-
mon functions you might find on a scientific calculator. This is especially true when using Jupyter
notebooks, since we don’t need to use the print() function to display the result of an expression
because the result of the last expression is always displayed automatically.

(2 + 2) * 3

12

# everything after a `#` is considered a comment and ignored by the␣


↪interpreter

# we can put multiple expressions in a cell


1 + 2 # result discarded
(continues on next page)

9
Python for Engineers

(continued from previous page)


3 + 4 # also discarded
5 + 6 # only the result of the last expression is displayed

11

# we can use print to manually show the result of any expression


print(1 + 2) # all of these are displayed
print(3 + 4)
print(5 + 6) # printing something and the last expression being␣
↪displayed are not quite the same

3
7
11

2 + * 3 # invalid input produces an error

Cell In[6], line 1


2 + * 3 # invalid input produces an error
^
SyntaxError: invalid syntax

Let’s calculate the area of a circle with radius 2.

2 * 2 * 3.14159

12.56636

1.3 Variables & Expressions

Without exponentiation, we have to write the radius two times, which is a bit annoying. To make
this easier, we can store the radius in a variable. This is done by using the assignment operator
=. While we are at it, we can also assign the value of 𝜋 in a variable, so we can reuse it later.
When we write multiple lines, the lines will be executed in order, and the result of the last line will
be displayed.

# assigns the value of the expression on the right to the variable on the␣
↪left
radius = 3 # ie `radius` is assigned the value 3, not the same as `=` in␣
↪math
pi = 3.14159
radius * radius * pi # use variables in place of values

28.27431

Now we can just change the value of the radius variable radius to 3 and the rest automatically
adapts to it.
Attempting to refer to an undefined variable is an error, variable names and all other identifiers are
case sensitive.

1.3. Variables & Expressions 10


Python for Engineers

x = undefined_variable # attempting to use an undefined variable␣


↪produces an error

------------------------------------------------------------------------
↪---
NameError Traceback (most recent call␣
↪last)
Cell In[9], line 1
----> 1 x = undefined_variable # attempting to use an undefined␣
↪variable produces an error

NameError: name 'undefined_variable' is not defined

x = Radius # we defined the radius with a lower case `r`, not an upper␣
↪case `R`

------------------------------------------------------------------------
↪---
NameError Traceback (most recent call␣
↪last)
Cell In[10], line 1
----> 1 x = Radius # we defined the radius with a lower case `r`, not␣
↪an upper case `R`

NameError: name 'Radius' is not defined

Variables can also be changed during runtime, they are called “variables” after all.

x = 2
print(x)
x = 3
print(x)

2
3

With 𝑟 and 𝜋 stored in variables, we can also calculate the circumference of the circle.

2 * radius * pi

18.849539999999998

We can also store the circumference in a variable. In general, everywhere you can write a value,
you can also write an expression. This is what makes programming so powerful.

circumference = 2 * radius * pi # and assignment has no result


circumference # so we write it out again (a single variable expression)␣
↪to display it

18.849539999999998

Checkup Question
Consider the following code snippet:

1.3. Variables & Expressions 11


Python for Engineers

radius = 5
Radius = 10
result = radius * 2
print(result)

What will be the output?


• A: 10
• B: 20
• C: An error because radius is defined twice.
• D: An error because Radius is not a valid variable name.
The correct answer is A: 10. While Radius is a valid variable name it does not refer to the same
variable as radius.

1.4 Types & Values

Unlike a normal calculator, we are not limited to numeric values. There are also other types of
values like text. In the programming context we generally call these strings because they are
strings of characters. Formally, we could say that a type is a set of values, i.e. the type integer is
the set of all whole numbers or string is the set of all sequences of characters. Informally the type
is the kind of value we are talking about just we might categorize other things.
Let’s use some strings to make the output a bit nicer. For different types there are different oper-
ations we can perform on them. We can do arithmetic with numbers, but we can’t do arithmetic
with strings. We can, however, concatenate them using the + (plus) operator. It is important to
note that even though we are using the same symbol, addition and concatenation are different
operations.

first_word = "Hello" # string literal values are enclosed in double␣


↪quotes

second_word = 'World' # or single quotes


combined_string = first_word + second_word
print(combined_string)

HelloWorld

# this does not work both operands need to be strings to perform␣


↪concatenation

print("circumference = " + circumference) # uncomment this line to see␣


↪the error

------------------------------------------------------------------------
↪---
TypeError Traceback (most recent call␣
↪last)

Cell In[15], line 2


1 # this does not work both operands need to be strings to␣
↪perform concatenation
----> 2 print("circumference = " + circumference) # uncomment this␣
(continues on next page)

1.4. Types & Values 12


Python for Engineers

(continued from previous page)


↪line to see the error

TypeError: can only concatenate str (not "float") to str

Just like in natural language, even if the syntax or grammar is correct, the sentence might still
not make sense. In Python something + something_else could be a syntactically correct
expression but if something is a number and something_else is a string it doesn’t make
sense and therefore we get a TypeError. To see what type a value has, we can use Python’s
built-in type() function.

# we can check the type of a value using the type() function


print(type(circumference)) # float is a type to store real numbers
print(type("circumference = ")) # this is a string

<class 'float'>
<class 'str'>

In some cases, we can convert values to another type, but this is not always without information
loss. For example, we know that we can convert a number to a string because that’s how we
represent numbers on paper. We perform this conversion by calling the type, in this case str, as
a function with the value we want to convert as its argument.

number_as_string = str(circumference)
print(number_as_string) # does the same as just printing the number␣
↪because printing needs to convert to string anyway
print(type(number_as_string)) # now it is a string

18.849539999999998
<class 'str'>

# conversion is not always possible


print(int("42")) # the string "42" can interpreted as an integer
print(int("hello")) # but "hello" does not make sense as an integer

42

------------------------------------------------------------------------
↪---
ValueError Traceback (most recent call␣
↪last)

Cell In[18], line 3


1 # conversion is not always possible
2 print(int("42")) # the string "42" can interpreted as an␣
↪integer
----> 3 print(int("hello")) # but "hello" does not make sense as an␣
↪integer

ValueError: invalid literal for int() with base 10: 'hello'

# if we convert the circumference to a string we can concatenate it to␣


↪the other string
print("circumference = " + str(circumference)) # now it works

1.4. Types & Values 13


Python for Engineers

circumference = 18.849539999999998

The four most important basic data types in Python are:


• int for integers meaning whole numbers
• float for floating point numbers an approximation of real numbers
• str for strings meaning sequences of characters
• bool for boolean or logical values (more on this later)
Every value or variable has a type. Since Python is a dynamically typed language, the type of a
variable can change and is determined automatically by the value it refers to. Even though types
might not be as visible and the rules surrounding them not as strict as in other languages, they
are still there.

print(type(2))
print(type(2.5))
print(type("Hello World"))
print(type(True))

<class 'int'>
<class 'float'>
<class 'str'>
<class 'bool'>

Checkup Question
What is the result of the expression "Radius = " + str(5.0)?
• A: An error: TypeError
• B: The string "Radius = 5.0"
• C: The string "Radius = str(5.0)"
• D: The number 5.0
The correct answer is B: The string "Radius = 5.0". Using the str() function on the number
5.0 converts it to a string which can be concatenated with the first string.

1.5 Functions

Just like we assigned values to variables to reuse them in expressions, we can define functions
to reuse entire expressions (and pretty much everything else). As a first intuition, we can think of
functions just like we know them from math. A function takes some input, substitutes it into an
expression and finally returns the expressions result.
We can define a function using the def keyword followed by its signature and definition. A func-
tions’ signature describes the amount and type of arguments a function takes, and what type of
value it returns.
For example the function 𝑓(𝑥) = 𝑥2 that maps a real number to another real numbers i.e. 𝑓 ∶ ℝ →
ℝ would look like this in Python:

1.5. Functions 14
Python for Engineers

# the definition consists of the def keyword, the function name and␣
↪parameter names

def f(x: float) -> float:


# the body of the function must be indented by 4 spaces (usually)
return x * x # return exits the function and _returns_ the provided␣
↪value

# everything not indented is outside the function


print(f(2.0))
print(f(4.2)) # and it returns the square of that value

4.0
17.64

While in mathematics it is very common to give functions single letter names, it is considered bad
practice in programming. A functions name should reflect what the function does. So in this case it
should not be named f but rather square. That way users of the function or someone reading the
code immediately can understand what the function does without having to look at the definition.
Functions in programming languages, also sometimes called procedures, can actually do much
more than mapping inputs to outputs. We will learn more about this in the next lecture.
Going back to our circle area example, we might define the function circle_area that computes
a circles’ area given its radius.

def circle_area(radius: float) -> float: # type annotations are optional␣


↪and can contradicted
# meaning we can call a function with an incompatible type with no␣
↪error
pi = 3.14159
return radius * radius * pi

print(circle_area(1.0))
print(circle_area(2.0))
# no error because python does not check types and `True` is interpreted␣
↪as the numeric value 1
print(circle_area(True))

3.14159
12.56636
3.14159

The type annotations are actually completely optional and not enforced, they are only there to help
humans understand whats going on. As far as Python is concerned we defined the function as
def circle_area(r): and the type annotations are just comments.
However, we do get errors when we pass to many or to few parameters to a function:

circle_area() # the function needs the radius to work

------------------------------------------------------------------------
↪---

TypeError Traceback (most recent call␣


↪last)
Cell In[23], line 1
(continues on next page)

1.5. Functions 15
Python for Engineers

(continued from previous page)


----> 1 circle_area() # the function needs the radius to work

TypeError: circle_area() missing 1 required positional argument: 'radius


↪'

circle_area(1, 2) # and it also does not know what to do with additional␣


↪parameters

------------------------------------------------------------------------
↪---
TypeError Traceback (most recent call␣
↪last)

Cell In[24], line 1


----> 1 circle_area(1, 2) # and it also does not know what to do with␣
↪additional parameters

TypeError: circle_area() takes 1 positional argument but 2 were given

Checkup Question
In the function definition def circle_area(radius: float) -> float:, what does ->
float indicate?
• A: The function requires exactly one input named float.
• B: It’s a comment ignored by Python.
• C: The function converts the input radius to a float.
• D: The function is expected to return a value of type float.
The correct answer is D: The function is expected to return a value of type float, since it is a
type hint, it must not necessarily return a float. In some sense B is also correct since it is an
annotation that Python does not really care about. However, unlike an actual comment the hinted
return type must be a valid expression.

1.6 Conditionals with if and else

The thing that makes programming much more powerful than a simple sequence of math equa-
tions, is the power to repeat or skip parts of the code based on values calculated in the program
itself.
On such “control structure” is the if - else construct which enables us to do something, some-
thing else or even nothing at all.
In this example, we want to print something different depending on the value of my_number.
The statements that are part of the if and else sections are grouped by indentation.
my_number = 2 # one equals sign is assignment, a statement (!)

if my_number == 2: # two equals signs is check for equality, a question␣


↪(?)
print("your number is 2") # this may or may not be executed

(continues on next page)

1.6. Conditionals with if and else 16


Python for Engineers

(continued from previous page)


if my_number > 3:
print("your number is greater than 3") # one or the other of these␣
↪two is executed

else:
print("your number is not greater than 3") # the other

your number is 2
your number is not greater than 3

In the context of functions, we can also see this as a way to define a function piecewise. For
example, we would define the absolute value of a number as follows:

𝑥 if 𝑥 ≥ 0
|𝑥| = abs(𝑥) = {
−𝑥 else

In Python it would look like this:

def absolute_value(x: float) -> float:


if x >= 0:
return x # we can have multiple return statements in a function
else:
return -x # but only one of them is ever executed

print(absolute_value(-1.0))
print(absolute_value(2.0)) # this is the absolute value function
print(abs(-3.0)) # python already has a built-in absolute value function

1.0
2.0
3.0

Checkup Question
What is true about expressions of the form [left side] = [right side]?
• A: Always evaluates to True or False.
• B: The left and right side must be equal.
• C: The left and right side will be equal.
• D: The left side must be a variable.
The correct answer is D: It is an assignment which is technically not an expression but a statement.
After it has been executed the variable on the left will have the value that the right side evaluated
to. If the right side depends on the variable the two sides will most likely not be equal afterward.

1.6. Conditionals with if and else 17


Python for Engineers

1.7 Operators

Python, like most programming languages, implements all the common arithmetical operators that
you already know, but there are also some more “computing” related operators that we will get to
in a bit.

1.7.1 Arithmetical Operators

Name Operator Example


Addition + 2 + 2 = 4
Subtraction - 5 - 3 = 2
Multiplication * 2 * 3 = 6
Division / 7 / 4 = 1.75
Integer Division // 7 // 4 = 1
Modulo % 7 % 4 = 3
Exponentiation ** 2 ** 3 = 8

For positive numbers the modulo operator is equivalent to the remainder of a division.

1.7.2 Relational Operators

Relational operators compare two numerical values and return what is called a boolean value. A
boolean value is always either True or False.

Name Operator Example


Equal == 2 == 2 -> True
Not Equal != 2 != 2 -> False
Greater > 2 > 3 -> False
Greater or Equal >= 3 >= 2 -> True
Less Than < 4 < 4 -> False
Less Than or Equal <= 4 <= 4 -> True

Note the difference between comparison == and assignment =, one checks for equality of the left
and right side, while the other assigns the right side to the left side.

1.7.3 Logical Operators

Lastly, there are the logical operators and, or and not. These are used to combine boolean
values. The effect of these can be seen in the so-called “truth tables”.

1.7. Operators 18
Python for Engineers

x y x and y x or y not x
False False False False True
False True False True True
True False False True False
True True True True False

print("will the student pass the course?")


attends_lecture = True
does_practice_problems = True
if attends_lecture and does_practice_problems:
print("good chance at passing!")
elif attends_lecture or does_practice_problems:
print("might have a shot")
else:
print("not looking good...")

will the student pass the course?


good chance at passing!

1.7.4 Operator Precedence

Just like we have rules like for the precedence of arithmetic operations like multiplication and
addition in mathematics, Python extends these rules to the new types of operators. Here is an
incomplete list from strongest to weakest precedence:
• () Parentheses
• ** Exponentiation
• *, /, //, % Multiplication and division related operators
• +, - Addition and subtraction
• <, <=, >, >=, ==, != Comparison operators
• not, and, or Logical operators
Checkup Question
What is the result of the expression 5 + 2 * 3 == 11?
• A: True
• B: False
• C: 21
• D: An error, because it is not a valid expression.
The correct answer is A: True. Just like in regular mathematics, multiplication is done before
addition. After that the comparison evaluates to True.

1.7. Operators 19
Python for Engineers

1.8 Examples

Now, to get a better understanding of how what you have learned can be used, let’s program
some simple examples.

1.8.1 Beam Deflection

A beam is rigidly supported at one end and a force 𝐹 is applied at the other free end like shown
in the figure below.

In you mechanics script you might find the following formula for this:
𝐹 2
𝑤(𝑥) = 𝑥 (3𝐿 − 𝑥)
6𝐸𝐼
Write a function to determine the deflection 𝑤 at position 𝑥 as a function of the applied force 𝐹
and length of the beam 𝐿. Assume that the stiffness 𝐸𝐼 = 1 is unity.

def beam_deflection(position: float, force: float, length: float) ->␣


↪float:
return force / 6 * position**2 * (3 * length - position)

print(beam_deflection(1, 3, 1))

1.0

1.8.2 Maximum

Write a function that takes two inputs and always returns the larger one.

def maximum(x: float, y: float) -> float:


if x > y:
return x
else:
return y

print(maximum(1, 2))
print(maximum(3, 4))
print(max(5, 6)) # again, python already has a built-in function for this

2
4
6

1.8. Examples 20
Python for Engineers

1.8.3 Price Calculation

A hardware store is selling bolts, they cost 20 cents a piece with bulk discounts when buying larger
quantities:
• 10% discount when buying 100 pieces or more
• 20% discount when buying 1000 pieces or more
• 30% discount when buying 10000 pieces or more
Write a function that calculates the price for a given quantity of bolts, taking into account the
discounts as they apply.

def calculate_price(count: int) -> float:


price_per_item = 0.2
discount = 0
if count >= 10000:
discount = 0.3
elif count >= 1000:
discount = 0.2
elif count >= 100:
discount = 0.1
return count * price_per_item * (1 - discount)

print(calculate_price(50))
print(calculate_price(200))
print(calculate_price(1000))
print(calculate_price(10000))

10.0
36.0
160.0
1400.0

1.9 Summary

In this lecture we learned about:


• Assigning variables and using them in expressions.
• Atomic data types: int, float, str, bool.
• Defining functions to reuse code somewhere else.
• Conditional code execution with if and else.
• Overview of mathematical, relational and logical operators.
• Practiced all of the above with some examples.

1.9. Summary 21
CHAPTER

TWO

MOVING ON

Open in JupyterLite: Solved Blank


In the last lecture we got started with programming in Python. In this lecture we will finish up the
basics with loops and more details about function.

2.1 Loops

In the last lecture we saw how if-else constructs are in some sense similar to piecewise func-
tions in mathematics. Like with the definition of the absolute value:

𝑥 if 𝑥 ≥ 0
abs(𝑥) = {
−𝑥 else

Which translates very closely to Python:

def absolute_value(x: float) -> float:


if x >= 0:
return x
else:
return -x
# in this special case we could also write it like this
# return x if x >= 0 else -x

2.1.1 for Loops

A similar comparison can be made with for loops and sums or products. Let’s say we want to
calculate the sum of the first 𝑛 square numbers. Mathematically would write this as:
𝑛
𝑓(𝑛) = ∑ 𝑖2
𝑖=1

For every 𝑖 from 1 to 𝑛 inclusive we want to add its square to our total sum. In Python this translates
to:

def sum_of_squares(n: int) -> int:


res = 0 # this is implicit in the summing operation, in python we␣
↪need to write it down
# everything in the line below has a correspondence to the␣
(continues on next page)

22
Python for Engineers

(continued from previous page)


↪mathematical formula
# we write `n + 1` because the interval made by range is exclusive /␣
↪open on the right

for i in range(1, n + 1):


# the `i**2` corresponds to the expression inside the sum, the␣
↪rest is implicit
res = res + i**2
# res += i**2 # equivalent to the line above
return res
# return sum(i**2 for i in range(1, n + 1)) # in this special case we␣
↪can also write it like this

print(sum_of_squares(1))
print(sum_of_squares(2))
print(sum_of_squares(3))

1
5
14

2.1.2 The range Function

We use the range function, to generate the series of values that we want to sum, or more generally,
iterate over. We get more information about a function by using Python’s built in help function:

help(range)

Help on class range in module builtins:

class range(object)
| range(stop) -> range object
| range(start, stop[, step]) -> range object
|
| Return an object that produces a sequence of integers from start␣
↪(inclusive)
| to stop (exclusive) by step. range(i, j) produces i, i+1, i+2, ...,
↪ j-1.
| start defaults to 0, and stop is omitted! range(4) produces 0, 1,␣
↪2, 3.
| These are exactly the valid indices for a list of 4 elements.
| When step is given, it specifies the increment (or decrement).
...

We have already seen that we can supply it with a pair of inclusive starting value and exclusive
stopping value:

# the stop value is exclusive so we need to add 1 to include 6


for i in range(1, 6 + 1):
print(i)

1
2
(continues on next page)

2.1. Loops 23
Python for Engineers

(continued from previous page)


3
4
5
6

If we only pass one argument to it, it will use that as the stopping value, implicitly set the starting
value to 0:

# starts at 0 and goes up to 4 for a total of 5 elements


for i in range(5):
print(i)

0
1
2
3
4

While it may be somewhat unintuitive, the fact that we start at zero but stop before 5 means that
we still end up with 5 numbers and 5 iterations.
Finally, instead of always incrementing by 1 from one number to the next we can also specify a
step size by passing a third value:

# count up in increments of 2
for i in range(0, 10, 2):
print(i)

0
2
4
6
8

If we only care about the fact that something is done a certain number of times, we can explicitly
ignore the iteration variable by assigning it to _:

# use an underscore to indicate that we want to ignore the value


# this is a convention and python does not care about it
for _ in range(3):
print("hello world")

hello world
hello world
hello world

2.1. Loops 24
Python for Engineers

2.1.3 Iterating Over a Collection

You might also remember that there is a special sigma notation for summing over a collection or a
set of values. Say we are solving a mechanics problem and we have a set of forces acting on an
object in the 𝑥 direction 𝐹𝑥 = {2 N, −3 N, 5 N}. We want to calculate the resulting force 𝐹res as:

𝐹res = ∑ 𝑓
𝑓∈𝐹𝑥

So we are summing over all elements 𝑓 in the collection or set of forces 𝐹𝑥 .


This notation also has a close equivalent in Python:

F_x = [2, -3, 5] # not actually a set, more on collections in the next␣
↪lecture

F_res = 0
for f in F_x:
F_res += f

print(F_res)
print(sum(F_x)) # we can also directly sum over collections

4
4

When we are writing for i in range(...) we are technically also looping over an iterable,
in this case the range iterable which gives us evenly spaced integers.
Checkup Question
Consider the following piece of code:

res = 0
for i in range(1, 5, 2):
res += i
print(res)

What will be the output?


• A: 4, B: 9, C: 10, D: 15
The correct answer is A: 4. Calling range(1, 5, 2) will generate the numbers from 1 to 5
(exclusive) in increments of 2. So the loop body will be executed for i = 1 and i = 3 which
sums to 4.

2.1.4 while Loops

Another kind of loop which is much more primitive and has no close equivalence in mathematics
is the while loop. As the name suggests, it loops “while” a condition is fulfilled, i.e. evaluates to
True. We can use this to emulate a for loop:

i = 0
while i < 5:
print(i)
(continues on next page)

2.1. Loops 25
Python for Engineers

(continued from previous page)


i += 1 # equivalent to `i = i + 1`
# don't forget to increment, otherwise the loop will run forever

0
1
2
3
4

This also nicely illustrates why you should prefer the for loop in this case: we can save ourselves
two lines of code and much more importantly two opportunities to make a mistake. Especially
forgetting to increment the counter and causing an infinite loop can be a very confusing bug.
Most of the time we want to loop over some collection or perform a known number of iterations,
in which case we use a for loop. There are much fewer cases where we use while loops, it is
generally in situations when it is not clear how many iterations something will take.
We might implement naive division by successive subtraction using while loops, if we didn’t know
about division we would have no way to know how often we need to subtract the divisor from the
dividend:

def naive_divide(number: int, divisor: int) -> int:


remainder = number
quotient = 0
# subtract the divisor from the remainder as long as it is larger␣
↪than the divisor

while remainder >= divisor:


remainder -= divisor
quotient += 1
return quotient

print(naive_divide(10, 3), 10 // 3)
print(naive_divide(7, 2), 7 // 2)
# print(naive_divide(7, -2), 7 // -2) # whoops, infinite loop

3 3
3 3

2.1.5 Executing Code Step by Step with a Debugger

As the flow of execution becomes more complex, it can quickly get hard to reason through the
code step by step in you head. In these cases, a debugger can be very helpful to do exactly that
while also showing us the current value of variables.
In VSCode, we can execute a cell in debug mode by pressing Ctrl + Shift + Alt + Enter
or by clicking the dropdown menu next to the run button and selecting “Debug Cell”. Executing
a cell in debug mode does not do anything by itself, to get control of the execution we must first
define a breakpoint where the interpreter will stop and wait for our instructions. Breakpoints can
be added by hovering over the blank space before a line and clicking on the red dot that appears.
The red dot signifies that a breakpoint has been set. Once execution has stopped we can use the
little menu that appears to either continue execution until the next breakpoint F5 or go to the next
line by pressing F10.

2.1. Loops 26
Python for Engineers

We will make use of the debugger in the next section to get a better understanding of what the
loop related keywords break and continue do.

2.1.6 break and continue

There are two additional keywords associated with loops: break and continue. break is used
to exit the loop prematurely, while continue is used to skip the rest of the current iteration and
continue with the next one. Let’s see how this works in practice:

# search for a number that is divisible by 7


for i in range(1, 10):
print(i) # set a breakpoint on this print statement, then start␣
↪debugging
# stop iterating after we found a number that is divisible by 7
if i % 7 == 0:
print("breaking:", i, "is divisible by 7")
break

1
2
3
4
5
6
7
breaking: 7 is divisible by 7

# print odd numbers from 1 to 10


for i in range(6):
# skip the current iteration if the number is even
if i % 2 == 0:
continue

print(i)

1
3
5

Checkup Question
Consider the following piece of code:

number = 5
while number > 0:
number = number - 1 # the same as `number -= 1`
if number % 3 == 0:
continue
print(number)
if number == 2:
break

What numbers will be in the output?


• A: 0, B: 1, C: 2, D: 3, E: 4, F: 5

2.1. Loops 27
Python for Engineers

The correct answers are C: 2 and E: 4. The loop counts down from the number 5. Before the first
printout the number is decremented to 4 skipping 5. The iteration for 3 is skipped with continue
because 3 is divisible by 3. After reaching and printing 2, the loop is exited with break.

2.2 More About Functions

We only got to know the very basics of functions in the last lecture, there is a lot more to know
about them. This is also not nearly everything, we will learn even more later on in the course.

2.2.1 Pure Functions

In the last lecture you already saw how functions in Python can behave just like functions in math
in that they are a deterministic mapping from a set of inputs to a set of outputs. But unlike functions
in math, this is not always true for functions in Python and other programming languages. Instead
of considering functions as mappings, we can think of them as a set of instructions or statements
which can potentially interact with the outside world and may ore may not return a value.
For example a functions return value could rely on a variable outside the function:

secret_variable = 42

# this functions return value depends on other things than the parameters,
↪ it depends on the value of `secret_variable`
def mystery_function() -> str:
return "the secret is " + str(secret_variable)

# just like with the definition when calling the function we also have to␣
↪use the brackets even without a parameter
print(mystery_function())
secret_variable = 5
print(mystery_function()) # same arguments, but different result

the secret is 42
the secret is 5

In this case, the functions return value can be different even if we call it with the same input.
The information can also flow in the other direction, meaning that instead of the outside world
effecting the functions result. The function can effect the outside world, for example by printing
something to the standard output.

def print_greeting() -> None: # this function takes no arguments, and␣


↪returns nothing
print("hello there")
# no `return` statement because we don't want to return anything
# the function ends automatically after the last line that is part of␣
↪the function

# calling the function a twice time makes a difference because the␣


↪function has "side effects"
print_greeting()
print_greeting() # produces a second hello there
(continues on next page)

2.2. More About Functions 28


Python for Engineers

(continued from previous page)

# this function has no side effects, calling it multiple times makes not␣
↪difference

x = sum_of_squares(3)
x = sum_of_squares(3)

hello there
hello there

Because of these deviations from the mathematical sense, functions in programming languages
are also sometimes called procedures to differentiate them from mathematical functions. If a pro-
cedure (or colloquially a function) behaves like a mathematical function, we call it a pure function.
You should prefer writing pure functions whenever possible because they are easier to reason
about since they are deterministic and (except for execution speed) it does not matter whether
you reuse the functions return value or call it another time since it cannot affect the outside world
and must always return the same value for the same inputs.

2.2.2 Docstrings

Another thing that helps other people understand your code, especially when reading it, are doc-
strings. They are written in triple quotes at the start of the function. They should briefly explain
what the purpose of a function is, what the input parameters are and what the function returns.
We will often omit these in the lecture in the name of brevity, but you should always use them when
writing anything serious.

def circle_area(radius: float) -> float:


# docstrings go at the beginning of the function
"""Computes the area of a circle with the given radius.

Args:
radius (float): The radius of the circle.

Returns:
float: The area of the circle.
"""
pi = 3.14159
return radius * radius * pi

Most code editors and IDEs will display this docstring when hovering over the function name. You
can also display them in the interpreter by using the help function:

help(circle_area) # displays signature and docstring

Help on function circle_area in module __main__:

circle_area(radius: float) -> float


Computes the area of a circle with the given radius.

Args:
radius (float): The radius of the circle.

(continues on next page)

2.2. More About Functions 29


Python for Engineers

(continued from previous page)


Returns:
float: The area of the circle.

2.2.3 Scopes and Variable Visibility

Every function has its own local “scope” or “namespace”. This means that variables defined inside
a function only persist during the function call and are not visible from the outside. The same is not
true in the other direction, however. A function can access variables defined outside the function
i.e. where it was defined. If a function defines a variable that is already defined in an outer scope,
the new variable will shadow the outer variable. This means that every reference to that variables
name inside the function will refer to inner variable and the outer variable will stay untouched.
Note: Unlike in many other programming languages, such as C or Java, if statements and loops
do not create a new scope. So unlike the braces that indentations replaces in Python, a new
level of indentation does not necessarily create a new scope.
Let’s take a look at some examples to get a better idea of how this works:
# try going through this step by step with a debugger and look at the␣
↪variables
outer_variable = 42
same_name = "hello"

def some_function() -> None:


print("outer_variable:", outer_variable) # accessing the outer␣
↪variable
local_variable = 5
same_name = "world" # this is a new variable that shadows the outer␣
↪variable
print("inside: same_name:", same_name) # this refers to the local␣
↪variable

# the function ends automatically after the last line that is part of␣
↪the function

some_function()
print("outside: same_name:", same_name) # this refers to the outer␣
↪variable, it is unchanged

print(local_variable) # this will raise an error because `local_


↪variable` is not defined in the outer scope

outer_variable: 42
inside: same_name: world
outside: same_name: hello

------------------------------------------------------------------------
↪---
NameError Traceback (most recent call␣
↪last)
Cell In[17], line 14
12 some_function()
13 print("outside: same_name:", same_name) # this refers to the␣
↪outer variable, it is unchanged
---> 14 print(local_variable) # this will raise an error because␣
(continues on next page)

2.2. More About Functions 30


Python for Engineers

(continued from previous page)


`local_variable` is not defined in the outer scope

NameError: name 'local_variable' is not defined

Checkup Question
Consider the following piece of code:

x = 5
def f():
x = 10
return x
x = 15
print(f(), x)

What will be the output?


• A: 15 15
• B: 10 15
• C: 10 10
• D: An error
The correct answer is B: 10 15. In the function x is set to 10 and returned which is printed. This
has no effect on the outer x. The outer x is simply set to 15 and printed afterward.

2.3 Examples

FizzBuzz is a simple programming exercise that is based on a children’s game. The rules are as
follows:
• Count up through the whole numbers
• If a number is divisible by 3 say Fizz instead of the number
• If a number is divisible by 5 say Buzz instead of the number
• If it is divisible by both 3 and 5 say FizzBuzz

for i in range(1, 30 + 1):


if i % 3 == 0 and i % 5 == 0:
print("FizzBuzz")
elif i % 3 == 0: # elif is short for else if and can be used for␣
↪other possible conditions besides the if
print("Fizz")
elif i % 5 == 0: # this means that only one of the if, else and else␣
↪ifs is executed
print("Buzz")
else:
print(i)

1
2
(continues on next page)

2.3. Examples 31
Python for Engineers

(continued from previous page)


Fizz
4
Buzz
Fizz
7
8
Fizz
Buzz
11
Fizz
13
14
FizzBuzz
16
17
Fizz
19
Buzz
Fizz
22
23
Fizz
Buzz
26
Fizz
28
29
FizzBuzz

2.3.1 Factorial

Write a function to implement the Factorial ! operator.


The factorial of a natural number is defined as follows:
𝑛
𝑛! ∶= ∏ 𝑖
𝑖=1

def factorial(n: int) -> int:


product = 1 # start with 1 as the product
for i in range(1, n + 1):
product = product * i
# product *= i
return product

print(factorial(0))
print(factorial(1))
print(factorial(3))
print(factorial(5))

1
1
6
120

2.3. Examples 32
Python for Engineers

Write a function that calculates the nth Fibonacci number.


The series of Fibonacci numbers is recursively defined as follows:

⎧ 0 𝑛=0
{
𝐹𝑛 = ⎨ 1 𝑛=1
{ 𝐹
⎩ 𝑛−1 + 𝐹𝑛−2 𝑛>1

# the basic approach


# 0 + 1 = 1
# 1 + 1 = 2
# 1 + 2 = 3
# 2 + 3 = 5
# 3 + 5 = 8

def fib(n: int) -> int:


f_before = 0 # these are local variables that only exist in the␣
↪function
f = 1
for _ in range(n):
# step through the iterations of the loop with a debugger to see␣
↪how the variables change
f_new = f + f_before
f_before = f
f = f_new
return f_before

print(fib(0))
print(fib(1))
print(fib(2))
print(fib(3))
print(fib(4))
print(fib(5))

0
1
1
2
3
5

# we can also implement it recursively, just like the original definition


def fib_rec(n: int) -> int:
if n == 0:
return 0
elif n == 1:
return 1
# you can put as many else ifs after the if for alternatives
# that are checked if all of the above are false
else:
return fib_rec(n - 1) + fib_rec(n - 2)

print(fib_rec(0))
print(fib_rec(1))
print(fib_rec(2))
print(fib_rec(3))
(continues on next page)

2.3. Examples 33
Python for Engineers

(continued from previous page)


print(fib_rec(4))
print(fib_rec(5))

0
1
1
2
3
5

2.3.2 Herons Algorithm

Heron’s algorithm is a very simple algorithm to get increasingly precise approximations for the
square root of a number.
A guess for the square root 𝑥𝑛 of the number 𝑎 is improved with the following computation:

𝑥𝑛 + 𝑎/𝑥𝑛
𝑥𝑛+1 =
2
Write a function to implement Heron’s algorithm, stop after the approximation is correct to at least
6 decimal places.

def heron_sqrt(num: float) -> float:


precision = 1e-6 # 1e-6 is short for 1*10^-6
root = 1 # start with a guess of 1
while abs(root - num / root) > precision: # keep improving the guess␣
↪until it is good enough

root = (root + num / root) / 2


return root

print(2**0.5)
print(heron_sqrt(2))

1.4142135623730951
1.4142135623746899

2.4 Summary

In this lecture we learned about:


• for loops and their similarities to sums and products in math.
• while loops and where you might need them.
• How to execute and reason through code step by step using a debugger.
• More details about functions like docstrings, scoping rules and the idea of pure functions.
• Practiced all of the above with some examples.

2.4. Summary 34
CHAPTER

THREE

BUILT-IN DATA STRUCTURES

Open in JupyterLite: Solved Blank


In this lecture we will get to know the four common data structures in Python. But what actually is a
data structure, and why would we need them? Previously, we already learned about atomic types
like int and float. They are atomic in the sense that we can not sensibly break them down into
smaller parts. A data structure also sometimes referred to as a compound type is a type that can
be broken down into smaller parts, eventually ending up with atomic types. Attempting to list the
reasons why we would want data structures already partially answers the question. Whenever
there is some structure to the data we are working with it should also be reflected in our code.
For example, let’s say we are writing a program to manage shopping lists. To find out how we
would model the shopping lists in our code we need to think about how a shopping lists is struc-
tured. One (of many) perspectives is that it is a list of pairs of the thing to buy and how many
of them to buy. We can represent the things to buy as strins and the amount as integers. For
the pair we could use another list of size 2, but a better idea is to use a tuple which functions
very similar but is conceptually different. In Python, we might then define the new shopping list
type as follows:

ShoppingList = list[tuple[str, int]]

An instance of a shopping list could then look like this:

shopping_list = [("milk", 2), ("apples", 5), ("bread", 1)]

But we are getting ahead of ourselves, let’s just start with lists, which is probably the most common
data structure.

3.1 Lists

If we want to store multiple elements in one place, we can use a list. Similar to lists in real life, a
list is an ordered collection of things. Here ordered means that each element has a clearly defined
position and two lists are only equal if they contain the same elements in the same order. This is
as opposed to sets which will introduce later. A list is created with square brackets.

list_1 = [1, 2, 3] # list of three numbers


# a list can contain multiple types of items and even other lists
list_2 = [42, "hello", "world", list_1, [4, 5, 6]]
empty_list = [] # empty list

(continues on next page)

35
Python for Engineers

(continued from previous page)


print(list_1) # the print function is able to print lists by itself
print(list_2)
print(type(list_2)) # the type of a list is list

print(len(list_1)) # get the length of a list with the len() function


print(len(list_2))

[1, 2, 3]
[42, 'hello', 'world', [1, 2, 3], [4, 5, 6]]
<class 'list'>
3
5

We can access the individual elements in a list using the [] operator and the index of the element
we want to access. Note that the indices start at 0, so if we want to get the second element of a
list, we have to use the index 1.

print(list_1[1])

This also works for strings, because strings are also kind of like lists of characters.
The square bracket operator generally works for all list-like things in Python. We will see more of
these later.

print("Hello World"[2])

It is also possible to index the elements from the end of the list by using negative indices.

# imagine list_1[-1] as a shorthand for list_1[len(list_1)-1]


print(list_1[-1]) # get the last element
print(list_1[-2]) # get the second to last element

3
2

Attempting to get an item at an index that does not exist will raise an IndexError.

list_1[3] # only indices 0, 1 and 2 are valid

------------------------------------------------------------------------
↪---
IndexError Traceback (most recent call␣
↪last)

Cell In[7], line 1


----> 1 list_1[3] # only indices 0, 1 and 2 are valid

IndexError: list index out of range

Just like strings, lists can be concatenated using the + operator

3.1. Lists 36
Python for Engineers

list_3 = [4, 5, 6]
list_4 = list_1 + list_3
print(list_4)
list_5 = list_1 * 5 # and multiplication as repeated addition also works
print(list_5)

[1, 2, 3, 4, 5, 6]
[1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3]

A lot of the times we need to do something with every element in a list.


Knowing the range() function and indexing using [] you might think of something like this:

course_participants = ["Alice", "Bob", "Charlie"]


# DON'T DO THIS
for i in range(len(course_participants)):
print("hello", course_participants[i])

hello Alice
hello Bob
hello Charlie

But there is a much cleaner way to do this, we can directly “iterate” through every item using a for
loop:

for participant in course_participants:


print("hello", participant)

hello Alice
hello Bob
hello Charlie

3.1.1 Operations on Lists

There are a few functions associated with lists that you should know about.

my_list = [2, 4, 6, 8]
my_list.append(1)
my_list.insert(1, 12) # new index of the element, and the element itself
print(my_list)

[2, 12, 4, 6, 8, 1]

my_list = [2, 4, 6, 8]

print(my_list.pop()) # remove last element


print(my_list)

print(my_list.pop(2)) # remove element with index


print(my_list)

3.1. Lists 37
Python for Engineers

8
[2, 4, 6]
6
[2, 4]

Be careful when inserting or removing elements at arbitrary positions in the list, this is a com-
paratively expensive operation as opposed to appending to or popping from the back. After an
insertion or deletion all elements after that position will need to be moved one position up or down
resulting in len(list)/2 move operations on average. Appending or popping from the back is
(for the most part) independent of the size of the list and therefore very fast.
chaos_list = [3, 4, 1, 2, 5]

chaos_list.sort() # sort list in place in ascending order


# chaos_list.sort(reverse=True) # descending order
print(chaos_list)
sorted_list = sorted([3, 4, 1, 2, 5]) # sorted returns a new list
print(sorted_list)

sorted_list.reverse() # reverse list in place


print(sorted_list)
print(list(reversed([1, 2, 3]))) # return new reversed list

[1, 2, 3, 4, 5]
[1, 2, 3, 4, 5]
[5, 4, 3, 2, 1]
[3, 2, 1]

In addition to single-element access, Python also offers the slicing operator [:] to easily get
sections of a list. The operator takes in two arguments: the start and end index of the slice. The
start index is inclusive, and the end index is exclusive. Both arguments are optional and default
to the start and end of the list, respectively.
other_list = ["I", "like", "lists", "in", "pthon", "!"]
print(other_list)

print(other_list[4]) # access single element with index


other_list[4] = "python" # fix the spelling error

# access a "slice" of the list, the starting index is inclusive, the␣


↪ending index is exclusive

print(other_list[2:5])

# if the start or end is not specified, its starts or end with the list
print(other_list[:3])
print(other_list[4:])
# this can also be used with negative indices, this leaves out the last␣
↪element of the list
print(other_list[:-1])

['I', 'like', 'lists', 'in', 'pthon', '!']


pthon
['lists', 'in', 'python']
['I', 'like', 'lists']
['python', '!']
(continues on next page)

3.1. Lists 38
Python for Engineers

(continued from previous page)


['I', 'like', 'lists', 'in', 'python']

Unlike regular indexing, slicing does not raise an IndexError if the indices are out of bounds.
Instead, the resulting slice only contains the elements that are in bounds.

other_list[3:100] # only return up to index 5

['in', 'python', '!']

There is also a variation on the slicing operator with two colons [::]. It lets us specify a step size
(just like with range()).

# print every second element starting at index 1 end ending a index 5


print(other_list[1:5:2])

# a negative step size can also be used to reverse a list


print(other_list[::-1])

['like', 'in']
['!', 'python', 'in', 'lists', 'like', 'I']

3.1.2 Multi-Dimensional Lists

Often we need to work with multidimensional lists, like with matrices or tables. Python does not
have multidimensional arrays, but we can just put lists in lists to make them ourselves.

list_2d = [
[1, 2, 3],
[4, 5, 6],
[7, 8, 9],
[10, 11, 12],
]
print(list_2d)

# indexing is done from the outer most list to the inner most list
print(list_2d[2]) # access third sublist
print(list_2d[2][1]) # access second element of third sublist
# this way the "coordinates" in the 2d list are addressed in reverse␣
↪order list_2d[y][x]

[[1, 2, 3], [4, 5, 6], [7, 8, 9], [10, 11, 12]]


[7, 8, 9]
8

Later we will see how NumPy makes working with such multidimensional data structures much
easier, especially mathematical operations.
Checkup Question
Consider the following code:

3.1. Lists 39
Python for Engineers

my_list = [1, 2, 3]
my_list.append(5)
print(my_list[3] + my_list[-3])

What will be the output?


• A: 4
• B: 6
• C: 7
• D: An error
The correct answer is C: 7. The first a[3] accesses the last element of the list, which is 5, and
the second a[-3] accesses the second element of the list, which is 2. So the result is 5 + 2 =
7.

3.1.3 List Examples

That was quite a lot of things, let’s see some of these in action to get an idea how they are used
in practice. We will write three functions: beam_deflection_list() which is supposed to
calculate the beam deflection at multiple points of a beam, is_sorted() which checks if a list
is sorted and fib_list() which generates the first N Fibonacci numbers.
# from last time ...
def beam_deflection(position: float, force: float, length: float) ->␣
↪float:
return force / 6 * position**2 * (3 * length - position)

def beam_deflection_list(positions: list[float], force: float, length:␣


↪float) -> list[float]:
res = []
for position in positions:
[Link](beam_deflection(position, force, length))
return res
# the above pattern is so common that there is a shorthand for it
# return [beam_deflection(position, force, length) for position in␣
↪positions] # list comprehension

print(beam_deflection_list([0, 1, 2], 6, 2))

[0.0, 5.0, 16.0]

def is_sorted(numbers: list[int]) -> bool:


# don't directly index into lists
# for i in range(len(numbers) - 1):
# if numbers[i] > numbers[i + 1]:
# return False
# even this is not an excuse for directly indexing into lists, we can␣
↪use zip
for left, right in zip(numbers[:-1], numbers[1:], strict=True): #␣
↪more about zip later
if left > right:
return False
(continues on next page)

3.1. Lists 40
Python for Engineers

(continued from previous page)


return True # if there is no unsorted pair, the list is sorted

print(is_sorted([1, 2, 3])) # True


print(is_sorted([1, 3, 2])) # False
print(is_sorted([])) # an empty list is always sorted

True
False
True

def fib_list(n: int) -> list[int]:


if n == 0:
return [0]

series = [0, 1]
for _ in range(n - 1): # use an underscore to make it clear that the␣
↪loop variable is not used
[Link](series[-1] + series[-2])
return series

print(fib_list(0))
print(fib_list(3))
print(fib_list(10))
print(fib_list(20))

[0]
[0, 1, 1, 2]
[0, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55]
[0, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1597,
↪ 2584, 4181, 6765]

3.2 Tuples

Tuples are in many ways very similar to lists, but despite their functional similarities, they are con-
ceptually very different. We create a tuple simply by separating the elements with a comma. Often
you will also see them enclosed in parentheses to prevent other syntax from taking precedence.

my_tuple = "hello", "world"


also_a_tuple = ("hello", "world") # these brackets are optional
print(type(my_tuple)) # tuple
print(type((1, 2))) # these parentheses are not optional because commas␣
↪also separate function arguments

<class 'tuple'>
<class 'tuple'>

Like lists, we can also read the elements with the [] operator. But unlike lists we can not write
them with the [] operator or append elements to them. This is because tuples are immutable,
they can not be changed after they are created.

3.2. Tuples 41
Python for Engineers

print(my_tuple[0]) # hello
my_tuple[0] = "goodbye" # this will raise an error, tuples are immutable
my_tuple.append("!") # so will this

hello

------------------------------------------------------------------------
↪---
TypeError Traceback (most recent call␣
↪last)

Cell In[22], line 2


1 print(my_tuple[0]) # hello
----> 2 my_tuple[0] = "goodbye" # this will raise an error, tuples are␣
↪immutable
3 my_tuple.append("!") # so will this

TypeError: 'tuple' object does not support item assignment

So why make this slightly awkward restricted version of a list, why not just use lists? The difference
lies in how you are supposed to use them: You should use tuples for possibly inhomogeneous
collections (multiple types) of fixed size, where it does not make sense to change one without
changing the other. And use lists for collections of arbitrary size made up of similar items that
may change.

# alice and her birth year belong together


person = ("Alice", 1999) # (name, birth year) different types

# it does not make sense to modify her birth year


# person[1] = 1998 # not possible

# it does however make sense to change the entire person however


person = ("Bob", 2002)

# list of rocks in my rock collection


rocks = ["smooth rock", "shiny rock", "big rock"] # same type

# it makes sense to change the rock a the second position


rocks[1] = "shinier rock"

# or add another item to the collection


[Link]("pointy rock")

A function taking a person tuple can reasonably expect, that it is always made up of the str
name and int birth year, but a function taking a list of rocks should expect any number of rocks
as strings. This is also reflected in the type annotation syntax.

tuple[str] # a tuple containing one string


tuple[str, int] # a tuple containing a string and an integer
list[str] # a list possibly containing any number of strings including␣
↪zero

list[str]

Tuples also implement an ordering so you can use <, >, <= and >= to compare and sort them.
The order is determined from the order of their elements, comparing the next pair of elements if

3.2. Tuples 42
Python for Engineers

the first ones are equal.

print((1, 3) > (1, 2)) # first pair of numbers tie, but 3 > 2 in the␣
↪second pair

print((1, 2) > (2, 1)) # first fair already decides 1 !> 2


print(sorted([(1, 3), (1, 2), (2, 1)])) # comparison enables sorting

True
False
[(1, 2), (1, 3), (2, 1)]

Checkup Question
Consider the following code:

my_tuple = (1, 2, 3)
my_tuple.pop(0)
print(my_tuple)

What will be the output?


• A: (1, 2, 3)
• B: (1, 2)
• C: (2, 3)
• D: An error
The correct answer is D: An error. Tuples do not have a pop() method, because they are im-
mutable. Attempting to call pop() on a tuple will raise an AttributeError.

3.3 Dictionaries

Dictionaries are a critical data structure that powers a lot of the underlying language. Similarly to
real world dictionaries, they map a discrete set of inputs which we call keys to associated outputs
which we call values. Most importantly, they do this very efficiently, so a dictionary is not just about
the association between the keys and values but also the fact that looking up a value by its key is
much faster than a naive implementation using a list of pairs.

en_to_de: dict[str, str] = { # english words (keys) to german words␣


↪(values)

"hello": "hallo",
"world": "welt",
}
print(type(en_to_de)) # dict
print(en_to_de["hello"], en_to_de["world"])

<class 'dict'>
hallo welt

We can check if a key is in a dictionary using the in operator, assigning to a key will create it if it
does not exist.

3.3. Dictionaries 43
Python for Engineers

from collections import Counter

sentence = "hello world!"


char_counts: dict[str, int] = {} # create an empty dictionary
for character in sentence:
# if character in char_counts:
# char_counts[character] += 1
# else:
# char_counts[character] = 1
char_counts[character] = char_counts.get(character, 0) + 1 # does␣
↪the same thing as the above
# get allows us to specify a default value if the key is not present

print(char_counts)
print(Counter("hello world!")) # the correct way to do it

{'h': 1, 'e': 1, 'l': 3, 'o': 2, ' ': 1, 'w': 1, 'r': 1, 'd': 1, '!': 1}
Counter({'l': 3, 'o': 2, 'h': 1, 'e': 1, ' ': 1, 'w': 1, 'r': 1, 'd': 1,
↪ '!': 1})

Iterating through a dictionary will yield its keys, we can iterate through the keys and values simul-
taneously using the items() method which will yield (key, value) tuples.

# DONT DO THIS
for char in char_counts:
print(f"Character '{char}' occurs {char_counts[char]} times in the␣
↪sentence")

print()
# DO THIS INSTEAD if you need both the key and the value
for char, count in char_counts.items():
print(f"Character '{char}' occurs {count} times in the sentence")

Character 'h' occurs 1 times in the sentence


Character 'e' occurs 1 times in the sentence
Character 'l' occurs 3 times in the sentence
Character 'o' occurs 2 times in the sentence
Character ' ' occurs 1 times in the sentence
Character 'w' occurs 1 times in the sentence
Character 'r' occurs 1 times in the sentence
Character 'd' occurs 1 times in the sentence
Character '!' occurs 1 times in the sentence

Character 'h' occurs 1 times in the sentence


Character 'e' occurs 1 times in the sentence
Character 'l' occurs 3 times in the sentence
Character 'o' occurs 2 times in the sentence
Character ' ' occurs 1 times in the sentence
Character 'w' occurs 1 times in the sentence
Character 'r' occurs 1 times in the sentence
Character 'd' occurs 1 times in the sentence
Character '!' occurs 1 times in the sentence

Here are some other commonly used dictionary methods:

3.3. Dictionaries 44
Python for Engineers

hotel_occupation: dict[str, int] = { # name, room number


"Alice": 101,
"Bob": 102,
"Charlie": 103,
}
print("guests:", hotel_occupation.keys())
print(f"{len(hotel_occupation)} rooms are occupied: {hotel_occupation.
↪values()}")

# Bob checks out


print("bob was in room", hotel_occupation.pop("Bob")) # pop removes a␣
↪key and returns its value
print(hotel_occupation)

# update a dictionary values from another dictionary


hotel_occupation.update({"Charlie": 202, "Eve": 203})
print(hotel_occupation)

guests: dict_keys(['Alice', 'Bob', 'Charlie'])


3 rooms are occupied: dict_values([101, 102, 103])
bob was in room 102
{'Alice': 101, 'Charlie': 103}
{'Alice': 101, 'Charlie': 202, 'Eve': 203}

Checkup Question
Consider the following code:

my_dict = {"a": 1, "b": 2}


my_dict["c"] = 3
my_dict.update({"a": 4, "d": 5})
print(my_dict.get("e", 0), my_dict["a"])

What will be the output?


• A: None 1
• B: 0 1
• C: 0 4
• D: An error
The correct answer is C: 0 4. The first part my_dict.get("e", 0) will return 0 because the key
"e" does not exist in the dictionary so the default value is used. The second part my_dict["a"]
will return 4 because we updated the value of the key "a" to 4.

3.3. Dictionaries 45
Python for Engineers

3.4 Sets

Sets are just like dictionaries (mostly) but without the values, they only have keys. They are pretty
much as mathematicians would define them, an unordered collection of unique elements.

my_set = {"A", "B", "C"}


print(type(set())) # we must create an empty set like this
print(type({})) # because this is an empty dictionary

print(my_set, len(my_set)) # notice how there is no order in a set

print("D" in my_set) # check if something is in a set


my_set.add("D") # add values using add
print("D" in my_set)
print(my_set)
my_set.add("D") # adding something already in the set does nothing
print(my_set)
my_set.remove("B") # take out values using remove
print(my_set)

<class 'set'>
<class 'dict'>
{'B', 'C', 'A'} 3
False
True
{'B', 'C', 'D', 'A'}
{'B', 'C', 'D', 'A'}
{'C', 'D', 'A'}

Unordered means that elements do not have a clearly defined position and two sets are considered
the same if they have the same elements as opposed to lists where this is not the case.

print({1, 2} == {2, 1}) # the same set


print([1, 2] == [2, 1]) # different lists

True
False

We can perform set operations on them:

set_a = {1, 2, 3}
set_b = {3, 4, 5}

print(set_a | set_b) # union


print(set_a & set_b) # intersection
print(set_a - set_b) # difference
print(set_a ^ set_b) # symmetric difference

{1, 2, 3, 4, 5}
{3}
{1, 2}
{1, 2, 4, 5}

Checkup Question
Consider the following code:

3.4. Sets 46
Python for Engineers

my_set = {1, 2, 3}
my_set.add(2)
my_set.remove(2)
print(my_set)

What will be the output?


• A: {1, 3}
• B: {1, 2, 3}
• C: {1, 3, 2}
• D: An error
The correct answer is A: {1, 3}. The first add(2) does not change the set because 2 is already
in it. The second remove(2) removes 2 from the set, leaving only 1 and 3.

3.5 More Examples

3.5.1 Why Sets over Lists?

Just like dictionaries, operations like test for membership and adding elements are much faster
than a naive implementation using lists.

from random import randint

even_numbers_list = list(range(0, 1000000, 2))


even_numbers_set = set(even_numbers_list)
random_numbers = [randint(0, 1000000) for _ in range(1000)]

%%time
# takes multiple seconds
even_count = 0
for num in random_numbers:
if num in even_numbers_list:
even_count += 1
print(even_count)

479
CPU times: user 6.02 s, sys: 14 ms, total: 6.04 s
Wall time: 5.64 s

%%time
# takes less than a tenth of a second
even_count = 0
for num in random_numbers:
if num in even_numbers_set:
even_count += 1
print(even_count)

3.5. More Examples 47


Python for Engineers

479
CPU times: user 661 μs, sys: 0 ns, total: 661 μs
Wall time: 602 μs

3.6 Dictionaries for Caching

Last lecture we saw that the recursive implementation for Fibonacci numbers was very slow and
limited to small inputs because it needed to recalculate the lower Fibonacci an exponentially in-
creasing number of times.

def naive_fib(n: int) -> int:


if n <= 1:
return n
return naive_fib(n - 1) + naive_fib(n - 2)

print(naive_fib(35)) # already takes a while

9227465

To get around this recomputing problem we can cache the result of previously computed Fibonacci
numbers in a dictionary. This way we can simply look up the previous two Fibonacci numbers
instead of recomputing them, thereby breaking the exponential growth of the number of recursive
calls.

fib_cache = {}

def cached_fib(n: int) -> int:


if n in fib_cache:
return fib_cache[n]
if n <= 1:
return n
result = cached_fib(n - 1) + cached_fib(n - 2)
fib_cache[n] = result
return result

print(cached_fib(100)) # just as fast as the iterative version (mostly)

354224848179261915075

This pattern, where we use cache the result of subproblems to speed up the calculation of larger
problems, has a name: it is called dynamic programming. If you ever actually need caching you
should not implement it yourself, instead you should use Python’s [Link] decorator
instead.

from functools import cache

@cache # this is a decorator, more on this later


def cached_fib(n: int) -> int:
if n <= 1:
return n
return cached_fib(n - 1) + cached_fib(n - 2)
(continues on next page)

3.6. Dictionaries for Caching 48


Python for Engineers

(continued from previous page)

print(cached_fib(100))

354224848179261915075

3.7 Summary

In this lecture we learned about Python’s four basic data structures: list, tuple, dict and set
and things we can do with them. As a general rule of thumb, here is when to use which:
• Use a list for similar items (type) of arbitrary size whose order matters, or checking for mem-
bership is not important.
• Use a tuple bundling related of different data types together, the size should usually be the
same from one instance to the next.
• Use dictionaries for efficient mapping from inputs to outputs of arbitrary size. If you have
a collection of dicts that all define the same keys you should probably use a namedtu-
ple or a dataclass instead, more on that in the lecture on classes and object-oriented
programming.
• Use a set to keep track of the existence of unique elements or to perform set operations.

3.7. Summary 49
CHAPTER

FOUR

NUMPY & MATPLOTLIB

Open in JupyterLite: Solved Blank


NumPy and Matplotlib are two of the most important libraries for scientific computing in Python.
What is a library? It is basically a collection of features that are not part of the language that we
can easily integrate into our own Python programs. This way, everyone doesn’t have to write
the same base functionality again before they can get to work. You should always prefer using
a library instead of implementing the functionality yourself. This saves a lot of time and reduces
room for errors, as libraries are usually thoroughly tested. It will also make your code easier to
understand for others that are familiar with the library.

4.1 Installing Packages

To make use of external libraries, we first need to install them. We can do this by invoking Pythons
package manager pip in your OS’s terminal (bash for Linux or PowerShell for Windows):

$ pip install numpy

You can also do this directly from the notebook by using a magic command:

%pip install numpy

...

This has the advantage that you do not need to make sure to install the package in the correct
environment since it will be installed in the same environment as the notebook is running in.
But beware where you install packages! In most cases you should not install packages for a
project into you main environment. This can quickly create version conflicts and might even break
some tools. Instead, you should create virtual environments for each project and install only the
packages you need for that project into that environment. We will talk about virtual environments
and how to manage them in a later lecture.

50
Python for Engineers

4.2 NumPy

As the name suggests, NumPy is a library for numerical computing in Python. We import it, using
an import statement. We can also import it under a short alias to avoid overly verbose code. Most
common libraries have an agreed upon alias that everyone uses. For NumPy, this is np. So we
import it like this:

# to import a library use the import statement


import numpy as np

# "as np" means we can write [Link] instead of [Link]

4.2.1 Note on Dot Notation

The . (dot) operator in Python and other programming languages is kind of similar to the “/” or
“" in your file systems. It is used to access members functions or variables of an object, similar
to how you would access a file in a folder. You can imagine NumPy as a folder and the array
function as a file in that folder. We can access the array function by typing [Link] (or
[Link] if you import NumPy as np). And similarly to how folders can also be nested in other
folders, objects can also be nested in other objects. NumPy, for example, contains a submodule
linalg which contains the function inv which can be accessed by [Link].
It is also possible to only import specific parts of a package using the from keyword. We will go
into more details about modules and importing later.

from [Link] import inv

If you have experience working in a command line environment you might know that you can use
something like dir or ls to list the contents of a folder. In Python, we can use the dir function
to do the same.

my_list = [1, 2, 3]
# dir(my_list) # here we can find the append method we talked about last␣
↪lecture

dir(my_list)[36:] # skip the first 36 special methods

['__subclasshook__',
'append',
'clear',
'copy',
'count',
'extend',
'index',
'insert',
'pop',
'remove',
'reverse',
'sort']

# we can also use it to see all numpy methods and submodules


print(len(dir(np))) # there is a lot of stuff in numpy
dir(np)[:10] # let's just look at the first 10

4.2. NumPy 51
Python for Engineers

536

['False_',
'ScalarType',
'True_',
'_CopyMode',
'_NoValue',
'__NUMPY_SETUP__',
'__all__',
'__array_api_version__',
'__array_namespace_info__',
'__builtins__']

4.2.2 NumPy Arrays

The bread and butter of NumPy is the NumPy array, it is different from Python’s list in that it
supports arithmetical operations and behaves more like a mathematical vector or matrix. We can
create an ndarray from a list by using the array function. While these arrays offer a lot of the
features that lists do, the two should not and can not be used interchangeably. Some notable
differences are that NumPy arrays are homogeneous, meaning that all elements must be of the
same type, their contents are stored in a contiguous block of memory (important for performance),
and they are fixed size, meaning we cannot append or remove elements from them.

# create a numpy array from a list


array = [Link]([1, 2, 3])
print(array)
print(type(array))
print([Link]) # all elements are of type int64
# [Link](2) # error: ndarrays do not have an append method

[1 2 3]
<class '[Link]'>
int64

As the name ndarray suggests, we can also create higher dimensional arrays. We do this by
passing in a nested list.

list_2d = [[1, 2, 3], [4, 5, 6], [7, 8, 9], [10, 11, 12]]

np_2d = [Link](list_2d)
print(np_2d)
print(np_2d.shape)
print(type(np_2d.shape), np_2d.shape[1]) # shape is tuple, access␣
↪elements with []

print("")

[[ 1 2 3]
[ 4 5 6]
[ 7 8 9]
[10 11 12]]
(4, 3)
<class 'tuple'> 3

4.2. NumPy 52
Python for Engineers

Unlike nested lists, the “sub-arrays” of a NumPy array must all have the same length.
[Link]([[1, 2], [3]]) # error: sub lists have different lengths

------------------------------------------------------------------------
↪---
ValueError Traceback (most recent call␣
↪last)
Cell In[8], line 1
----> 1 [Link]([[1, 2], [3]]) # error: sub lists have different␣
↪lengths

ValueError: setting an array element with a sequence. The requested␣


↪array has an inhomogeneous shape after 1 dimensions. The detected␣
↪shape was (2,) + inhomogeneous part.

These homogeneity restrictions (while necessary for the performance benefits) also benefits us in
other ways like allowing more flexible slicing.
print(np_2d) # cleaner formatting
print(np_2d[2, 1]) # easier access [y, x] instead of [y][x]
print(np_2d[1:3, 1:3]) # access a rectangular subsection
print(np_2d[1]) # access second row
print(np_2d[:, 1]) # access second column

[[ 1 2 3]
[ 4 5 6]
[ 7 8 9]
[10 11 12]]
8
[[5 6]
[8 9]]
[4 5 6]
[ 2 5 8 11]

list_3d = [[[1, 2], [3, 4]], [[5, 6], [7, 8]]]

np_3d = [Link](list_3d) # we can have any number of dimensions


print(np_3d.shape)

(2, 2, 2)

4.2.3 Vectorized Computations

A lot of the time when working with numerical data, we want to apply the same operation to all
elements of an array or perform element-wise operations between two arrays. The homogeneous
nature of ndarrays not only makes this very performant, it also makes it very easy to read and
write. Most arithmetic operations that we can perform on single numbers (or pairs of numbers)
we can also perform on arrays without explicitly looping over them.
# scalar operations on an array
array2 = array * 2 + 2
print(array2)
(continues on next page)

4.2. NumPy 53
Python for Engineers

(continued from previous page)

# element wise operations


print(array * array2)

[4 6 8]
[ 4 12 24]

When doing element wise operations, the shapes must be compatible. This usually means that
they must either have the same length or one of them must be a scalar in which case scalar will be
“broadcasted” to the shape of the other array. Understanding these broadcasting rules for higher
dimensional arrays can be a bit confusing at first, but can make certain operations very convenient
to express.

[Link]([1, 2, 3]) + [Link]([4, 5]) # incompatible shapes

------------------------------------------------------------------------
↪---
ValueError Traceback (most recent call␣
↪last)
Cell In[12], line 1
----> 1 [Link]([1, 2, 3]) + [Link]([4, 5]) # incompatible shapes

ValueError: operands could not be broadcast together with shapes (3,)␣


↪(2,)

Theses implicit vectorized operations mean that all function that only use arithmetic operations
naturally work on arrays as well.

# implicitly also (a: float | [Link], b: float | [Link]) -> float␣


↪| [Link]
def add_and_square(a: float, b: float) -> float:
# could also be (a: [Link], b: [Link]) -> [Link]
return (a + b) ** 2

print(add_and_square(2, 5))
print(add_and_square([Link]([1, 2, 3]), [Link]([4, 5, 6]))) # also␣
↪works with arrays

49
[25 49 81]

NumPy also comes with all the usual mathematical constants and functions, like the square root,
power, and trigonometric functions. If they are used on a NumPy array (or any other iterable), the
function is applied to each element.

vec = [Link]([1, 2, 3, 4])

squares = vec**2
print(squares)
print([Link](squares))
print([Link]([0, [Link] / 2]))

4.2. NumPy 54
Python for Engineers

[ 1 4 9 16]
[1. 2. 3. 4.]
[0. 1.]

4.2.4 Array Creation Functions

There are a few functions that create special arrays that are commonly needed for doing vec-
torized computations. The most common ones are linspace, arange, zeros, ones, and
[Link].

# linspace creates a series of evenly spaced values


print([Link](0, 10, 41)) # start, stop, count
print([Link](0, 24.5)) # the default count is 50

[ 0. 0.25 0.5 0.75 1. 1.25 1.5 1.75 2. 2.25 2.5 2.75


3. 3.25 3.5 3.75 4. 4.25 4.5 4.75 5. 5.25 5.5 5.75
6. 6.25 6.5 6.75 7. 7.25 7.5 7.75 8. 8.25 8.5 8.75
9. 9.25 9.5 9.75 10. ]
[ 0. 0.5 1. 1.5 2. 2.5 3. 3.5 4. 4.5 5. 5.5 6. 6.5
7. 7.5 8. 8.5 9. 9.5 10. 10.5 11. 11.5 12. 12.5 13. 13.5
14. 14.5 15. 15.5 16. 16.5 17. 17.5 18. 18.5 19. 19.5 20. 20.5
21. 21.5 22. 22.5 23. 23.5 24. 24.5]

# arange is similar but instead of specifying how many values are in the␣
↪array, the distance between them is specified
print([Link](1, 10, 1)) # almost equal to [Link](range(1,10,1))
print([Link](1, 10, 0.5))

[1 2 3 4 5 6 7 8 9]
[1. 1.5 2. 2.5 3. 3.5 4. 4.5 5. 5.5 6. 6.5 7. 7.5 8. 8.5 9. 9.
↪5]

# zeros and ones create arrays of the specified size made up of 0s and 1s␣
↪respectively
print([Link](5))
print([Link](10))
# the shape can also be a tuple to create multi-dimensional arrays
print([Link]((2, 3)))

# random creates arrays of specified size with a random value from a␣


↪uniform distribution between 0 and 1
print([Link](8))
print([Link](4) * 10) # generate 4 random numbers between 0␣
↪and 10

[0. 0. 0. 0. 0.]
[1. 1. 1. 1. 1. 1. 1. 1. 1. 1.]
[[1. 1. 1.]
[1. 1. 1.]]
[0.38528082 0.58261164 0.24707753 0.41324455 0.89771061 0.35981933
0.78797988 0.48444789]
[2.45060753 4.63820712 0.14009443 0.69574187]

4.2. NumPy 55
Python for Engineers

Checkup Question
Consider the following code:

a = [Link]([1, 2])
b = [Link]([3, 4])
print(a + b)

What will be the output?


• A: [4, 6]
• B: [1, 2, 3, 4]
• C: [[1, 2], [3, 4]]
• D: An error
The correct answer is A: [4, 6]. Using the + operator on two NumPy arrays will do an element-
wise addition.

4.2.5 Linear Algebra

NumPy also implements a lot of linear algebra operations. There is no special type for matrices or
vectors in NumPy, it is the operation that defines whether an ndarray is interpreted as a matrix
/ vector or something list-like.
To multiply two matrices or vectors together, we could use the [Link](A, B) function. But
just like we can use ** instead of pow() we can also use the @ operator instead of [Link]()
which makes the code a lot more readable. When doing a matrix vector operation, NumPy will
automatically interpret the 1D array as a 2D array with one row or column, depending on the order
of the operands. This also means that multiplying two vectors, using the matrix multiplication
operator, will do the same as the dot product.

A = [Link]([[1, 2], [-1, 4]])


print("A=")
print(A)
B = [Link]([5, 6])
print("B=")
print(B)
print("A*B=")
print(A @ B) # use the @ symbol for matrix multiplication
print("B*A=")
print(B @ A) # numpy automatically transposes row vectors for␣
↪multiplication

A=
[[ 1 2]
[-1 4]]
B=
[5 6]
A*B=
[17 19]
B*A=
[-1 34]

4.2. NumPy 56
Python for Engineers

F = [Link]([[1, 1], [1, 0]]) # fibonacci matrix


# notice how the top right corner of F^n is the nth fibonacci number
print("F=")
print(F)
print("F*F=")
print(F @ F)
print("F^3=")
print(F @ F @ F)
print("F^5=")
print([Link].matrix_power(F, 5))

F=
[[1 1]
[1 0]]
F*F=
[[2 1]
[1 1]]
F^3=
[[3 2]
[2 1]]
F^5=
[[8 5]
[5 3]]

And of course we can also easily invert matrices and solve linear systems of equations.
# solve a linear system of equations in the form Ax = B
print([Link](A) @ B) # using inverse
print([Link](A, B)) # using solve

[1.33333333 1.83333333]
[1.33333333 1.83333333]

You should generally prefer [Link] over multiplying the inverse, as it is more nu-
merically stable and faster.
NumPy has functions for most other linear algebra operations you know about, like the determi-
nant, inverse, transpose, eigenvalues / eigenvectors or the singular value decomposition.
print([Link](A)) # determinant of a
print([Link](A)) # inverse of a
print([Link](A))
print(A.T) # if you matrix is already a numpy array should use .T to get␣
↪a transposed view

6.0
[[ 0.66666667 -0.33333333]
[ 0.16666667 0.16666667]]
[[ 1 -1]
[ 2 4]]
[[ 1 -1]
[ 2 4]]

# eigenvalues and vectors


print([Link](A))
(continues on next page)

4.2. NumPy 57
Python for Engineers

(continued from previous page)


print([Link]([Link]([3, 5]))) # diag creates a diagonal␣
↪matrix
print([Link](A)) # eigenvalues and vectors at the same time in a␣
↪namedtuple
vals, vecs = [Link](A) # destructure the tuple into the values␣
↪and vectors
print(vecs[:, 0], vecs[:, 0]) # the eigenvectors are columns of the␣
↪matrix

[2. 3.]
[3. 5.]
EigResult(eigenvalues=array([2., 3.]), eigenvectors=array([[-0.89442719,
↪ -0.70710678],
[-0.4472136 , -0.70710678]]))
[-0.89442719 -0.4472136 ] [-0.89442719 -0.4472136 ]

# singular value decomposition


print([Link](A))

SVDResult(U=array([[ 0.41785681, 0.9085129 ],


[ 0.9085129 , -0.41785681]]), S=array([4.49661478, 1.33433712]),␣
↪Vh=array([[-0.10911677, 0.99402894],
[ 0.99402894, 0.10911677]]))

Checkup Question
Consider the following piece of code:

A = [Link]([[1, 2], [3, 4]])


B = [Link]([[5, 6], [7, 8]])
print((A * B)[0, 0])

What will be the output?


• A: 5
• B: 19
• C: A ShapeError
• D: An IndexError
The correct answer is A: 5. The * operator does an element-wise multiplication of the two arrays.
The upper left element is then retrieved using the [0, 0] index. Matrix multiplication is performed
using the @ operator.

4.2. NumPy 58
Python for Engineers

4.2.6 Examples

Let’s look at a few examples to see how we might use the features we just learned about in prac-
tice.

Vectorized Square Numbers

Write a function that calculates the first 𝑛 square numbers and returns them as a NumPy array
without using a loop.

def square_numbers(n: int) -> [Link]:


return [Link](1, n + 1) ** 2

print(square_numbers(10))

[ 1 4 9 16 25 36 49 64 81 100]

Beam Deflection NumPy

In the first lecture, we wrote a function to calculate the beam deflection at a given point for a given
load and beam length:

def beam_deflection(position: float, force: float, length: float) ->␣


↪float:
return force / 6 * position**2 * (3 * length - position)

Exploit NumPy’s vectorized operations to compute the deflection at 10 equally spaced points
along a beam of length = 5 with a load of force = 10.

length = 5
positions = [Link](0, length, 10)
deflections = beam_deflection(positions, 10, length)
deflections

array([ 0. , 7.43026978, 28.57796068, 61.72839506,


105.16689529, 157.17878372, 216.04938272, 280.06401463,
347.50800183, 416.66666667])

Matrix Fibonacci Numbers

Write a function that calculates the 𝑛-th Fibonacci number repeated multiplication of the Fibonacci
matrix F:
1 1
F=[ ]
1 0

def fib_mat(n: int) -> int:


F = [Link]([[1, 1], [1, 0]])

(continues on next page)

4.2. NumPy 59
Python for Engineers

(continued from previous page)


res = [Link](2, dtype=np.int32) # start with identity matrix
for _ in range(n):
res @= F # equivalent to `res = res @ F`

return res[0, 1]
# or using matrix power
# return [Link].matrix_power(F, n)[0, 1]

print(fib_mat(0))
print(fib_mat(1))
print(fib_mat(2))
print(fib_mat(3))
print(fib_mat(25))

0
1
1
2
75025

Statics

Find the rope tension forces 𝑆𝐿 and 𝑆𝑅 in the arrangement below:

We can derive two equations by considering the force equilibrium conditions at the center node:
!
Σ𝐹𝑥 = sin(−45°) ⋅ 𝑆𝐿 + sin(60°) ⋅ 𝑆𝑅 = 0
!
Σ𝐹𝑦 = cos(−45°) ⋅ 𝑆𝐿 + cos(60°) ⋅ 𝑆𝑅 − 𝐹 = 0

Assuming that 𝐹 = 𝑚𝑔 = 10 reformulate the equation in matrix form and determine 𝑆𝐿 and 𝑆𝑅
using NumPy.

A = [Link]([
[[Link]([Link](-45)), [Link]([Link](60))],
[[Link]([Link](-45)), [Link]([Link](60))],
])
b = [Link]([0, 10])
S_L, S_R = [Link](A, b)
print(f"{S_L = :.3f}, {S_R = :.3f}")

4.2. NumPy 60
Python for Engineers

S_L = 8.966, S_R = 7.321

Eigendecomposition Fibonacci Numbers

Finally, let’s try to use some more in-depth mathematics knowledge to solve the problem with
fewer matrix multiplication operations than our simple solution. Some computations, like expo-
nentiation, can be performed more efficiently in a basis of eigenvectors.

M𝑛 = V ⋅ D𝑛 ⋅ V−1

Where D is a diagonal matrix with the eigenvalues of M and V a matrix whose columns are made
up of the respective eigenvectors of M.
Remember that the powers of a diagonal matrix can be computed by simply computing the powers
of the diagonal elements.

fib_mat = [Link]([[1, 1], [1, 0]])


eigenvalues, eigenvectors = [Link](fib_mat)
# notice how the eigenvalues of the fibonacci matrix are the two golden␣
↪ratios
print(eigenvalues)
# the eigenvectors are the columns of this matrix, this is already our V
print(eigenvectors)

V = eigenvectors
V_inv = [Link](V)
D = [Link](eigenvalues)

# this is using the same math as the exercise from the previous lecture
def fib_eig(n: int) -> int:
# round to account for floating point inaccuracies
return round((V @ [Link](eigenvalues**n) @ V_inv)[0, 1])

print(fib_eig(0))
print(fib_eig(1))
print(fib_eig(2))
print(fib_eig(3))
print(fib_eig(25))

[ 1.61803399 -0.61803399]
[[ 0.85065081 -0.52573111]
[ 0.52573111 0.85065081]]
0
1
1
2
75025

4.2. NumPy 61
Python for Engineers

4.3 Matplotlib

Just like with NumPy (and many other libraries), the primary use of the library is in the name. In
the case of Matplotlib it is creating all kinds of plots.

import [Link] as plt

[Link]("[Link]")

4.3.1 Basic Plotting

Plotting in Matplotlib is very similar to how you would plot with pen and paper. Simply give the plot
function a series of x and y values, and it will draw and connect them with lines.

np.set_printoptions(precision=4, threshold=20) # nicer printout in static␣


↪version

# x = [Link]([-2,-1,0,1,2]) # the segments are clearly visible when␣


↪only plotting a few points
xs = [Link](-2, 2)
ys = xs**2 # this will square each element of xs
print("xs:", xs)
print("ys:", ys)
[Link](xs, ys);

xs: [-2. -1.9184 -1.8367 ... 1.8367 1.9184 2. ]


ys: [4. 3.6801 3.3736 ... 3.3736 3.6801 4. ]

There are a lot of ways to customize plots, here are a few examples. Check out the cheat sheet
for an overview of the most important kind plots available and how to customize them.

fig, ax = [Link]()

(continues on next page)

4.3. Matplotlib 62
Python for Engineers

(continued from previous page)


xs = [Link](0, 2 * [Link])

[Link](xs, [Link](xs), label="sin(x)")


# if no color is specified pyplot will automatically cycle through a list␣
↪of colors
[Link](xs, [Link](xs), label="cos(x)")
[Link](xs, [Link](xs) / 500, color="purple", linestyle="dashdot", label=
↪"exp(x)/500")
# similar to [Link] but just uses unconnected markers
[Link]([0, [Link], 2 * [Link]], [0, 0, 0], color="red", marker="x",␣
↪label="zeros of sin(x)")
[Link]()

ax.set_xlabel("x [s]")
ax.set_ylabel("y [m]")
ax.set_title("Example of Matplotlib")
[Link]("[Link]") # save the figure to a {png, svg, pdf, ...}␣
↪file
# [Link](["a","b","c","d"]) # override the labels set in the plot␣
↪commands

Matplotlib also supports a lot of other types of plots, like bar charts, box plots, vector plots, and
many more. If you can think of it, it’s probably already implemented. If you need to use them, you
will find lots of examples and help in the documentation online.

4.3. Matplotlib 63
Python for Engineers

4.3.2 Subplots

Sometimes we want to plot multiple things in one figure instead of having them each in their own
figure. For this, we can use [Link](), the function returns the “handles” for the figure
and the different plots in a list. Instead of calling the plot functions with plt, we call them on the
handles of the different plots called axes, you can also work this way if you only have a single
figure.
When working with axes instead of plt, some commands change slightly:
[Link]() -> ax.set_title()
[Link]() -> ax.set_xlabel()
[Link]() -> ax.set_xlim()
... -> ...
If you are unsure, just Google it, look at the cheat sheet or search the Matplotlib documentation.

fig, axes = [Link](2, 2)

xs = [Link](0, 1) # some dummy data


ys = xs**2

axes[0, 0].plot(xs, ys) # top left corner


axes[0, 1].plot(-xs, ys) # top right corner
axes[1, 0].plot(xs, -ys) # bottom right corner
axes[1, 1].plot(-xs, -ys) # bottom right corner

axes[0, 0].set_xlabel("xlabel")
axes[1, 0].set_title("axes title")
[Link]("figure title");

4.3. Matplotlib 64
Python for Engineers

4.3.3 Examples

It is time for more examples. Plot the following functions in the range [−3, 3] and properly annotate
the plot:
1
𝑓1 (𝑥) = 2𝑥 − 𝑥3
4
𝑥
𝑓2 (𝑥) = sin(𝑒 )
1
𝑓3 (𝑥) =
𝑥
Make sure the plot looks correct and the region shown is sensible.
xs = [Link](-3, 3, 200)

f1s = 2 * xs - 0.25 * xs**3


[Link](xs, f1s, label="$f_1(x)$") # we can use latex in labels and␣
↪other matplotlib text

f2s = [Link]([Link](xs))
[Link](xs, f2s, label="$f_2(x)$")

f3s = 1 / xs
f3s[len(xs) // 2] = [Link]
# a nan in the y data will lead to a break in the lines
[Link](xs, f3s, label="$f_3(x)$")
# you can also produce the nan by using an odd number of points which␣
↪will contain a 0
# the division by zero will then also lead to a nan in the correct place

[Link]()
[Link](-5, 5) # we can also override the automatic range
[Link]("$x$")
[Link]("$f_n(x)$");

And for something slightly less “made-up”, plot the deflection of the beam from the previous sec-
tion along its length and annotate the plot. The force is applied in the negative y-direction, so the

4.3. Matplotlib 65
Python for Engineers

beam should deflect downwards.

positions = [Link](0, length)


[Link](positions, -beam_deflection(positions, 10, length))
[Link]("position")
[Link]("deflection");

Checkup Question
Consider the following code:

xs = [Link](-2, 2, 8)
[Link](xs, xs**2)
[Link](xs, 1/xs)

What is problematic about the resulting plot?


• A: The discontinuity is handled improperly
• B: The plot is not annotated properly
• C: The plotting resolution is too low
• D: The y-ranges of the functions are too different
The correct answers are A, B, and C. Matplotlib will connect the two points on either side of the
discontinuity leading to a misleading plot, there is also not annotations like axis labels legend or title
to explain what is being plotted and finally, the linspace with only 8 points means the resolution
is so low that the individual line segments are visible. The last point is not a problem, the parabola
is in the 𝑦-interval [0, 4] while the largest magnitude point of the hyperbola is ∼ 1/0.25 = 4 which
is on the same scale.

4.3. Matplotlib 66
Python for Engineers

4.4 Summary

In this lecture we learned about:


• How to install packages with pip.
• The NumPy library and common operations with it’s ndarray type.
• How use NumPy for linear algebra operations.
• The basics of plotting with Matplotlib.

4.4. Summary 67
CHAPTER

FIVE

STRINGS, FILES & FILE SYSTEMS

Open in JupyterLite: Solved Blank


In this lecture we are going to learn a bit about working with strings, files and file systems in Python.
To do this we will write an example project that requires us to use all of these features.

5.1 Task Description

As a diligent student, you have always taken notes in every class, and because you really like the
Markdown format, you have made all your notes in Markdown (the same markup language that is
used for these Jupyter notebooks). Unfortunately you are not that good at keeping things in order,
and you always have trouble finding the files where your took notes. So to get things in order, we
are going to write a program that will search for all markdown files, look for all headings in these
files, and then make a table of contents file that lists all the headings and links to them. To make it
a bit more interesting, we only want to list top and second level headings in the table of contents.
Here is an example to better illustrate what we want to achieve:
This is what the file system looks like:

notes
├── [Link]
│ > # Important Stuff
│ > ## How to Sign up for Exams
├── math
│ └── [Link]
│ > # My Math Notes
│ > not part of a heading
│ > #also not part of a heading
│ >
│ > ## Linear Algebra
│ > ### Matrices
│ > ## Differential Equations
└── mechanics
└── mechanics-1
├── [Link]
│ > # Mechanics 1 Exercises
└── [Link]
> # Mechanics 1
> ## Statics

And this is what the table of contents file should look like:

68
Python for Engineers

# Table of Contents

- notes/
- [[Link]](<notes/[Link]>)
- [Important Stuff](<notes/[Link]#important-stuff>)
- [How to Sign up for Exams](<notes/[Link]#how-to-sign-up-
↪for-exams>)
- math/
- [[Link]](<notes/math/[Link]>)
- [My Math Notes](<notes/math/[Link]#my-math-notes>)
- [Linear Algebra](<notes/math/[Link]#linear-algebra>)
- [Differential Equations](<notes/math/[Link]#differential-
↪equations>)
- mechanics/
- mechanics 1/
- [[Link]](<notes/mechanics/mechanics 1/[Link]>)
- [Exercises](<notes/mechanics/mechanics 1/[Link]#exercises>
↪)
- [mechanics [Link]](<notes/mechanics/mechanics 1/mechanics [Link]>)
- [Mechanics 1](<notes/mechanics/mechanics 1/mechanics [Link]
↪#mechanics-1>)
- [Statics](<notes/mechanics/mechanics 1/mechanics [Link]#statics>
↪)

To be able to complete this task, we will need to learn a few things:


• How to work with the file system
– How to represent a path to a file or directory
– How to figure out what type a file is
– How to traverse and search through subdirectories
• How to read and write files
• How to work with the file contents
– How to parse the contents of the file
– How to programmatically create the output file

5.2 Working with the File System

The file system on most modern computers is what we could call a hierarchical file system. And as
the name suggests this means there is some hierarchy to everything. More precisely, this means
that the files don’t just lay around in some big metaphorical pile of files, but they are structure in
directories and subdirectories. To uniquely refer to a file, we specify all the steps we need to take
to get to it, we call these steps the “path” to the file. Take this path for example:
C:\Users\Alice\Documents\[Link]
It tells us that we need to start at the C drive, then go to the Users directory, the Alice directory,
the Documents directory where we will find a file called [Link]. The individual steps of the
path are separated by a backslash \ so if we were given this path as a string and wanted to get
the individual components we could simply split it at the backslashes:

5.2. Working with the File System 69


Python for Engineers

path_string = r"C:\Users\Alice\Documents\[Link]" # use a raw string to␣


↪prevent escape characters

path_string.split("\\") # "\\" is a string containing a single backslash

['C:', 'Users', 'Alice', 'Documents', '[Link]']

But now when we go to a different operating system like Linux, the equivalent path would look
something like this:
/home/alice/documents/[Link]
There is a different path separator and the path starts at the root directory / instead of a drive
letter like on Windows. So previous example code would not work on Linux. Because of this and
many other reasons you should not directly work with paths as string and manually manipulate
them. To write robust platform independent code, we can use the pathlib module that is part of
Python’s standard library. It provides a class called Path that implements many useful methods
related to working with paths.

from pathlib import Path, PureWindowsPath

# todo_path = Path(r"C:\Users\Alice\Documents\[Link]") # does not work␣


↪properly on linux
# you wouldn't normally use PureWindowsPath, this is just so this little␣
↪demo also works on linux
todo_path = PureWindowsPath(r"C:\Users\Alice\Documents\[Link]") #␣
↪create a path from a string
# path = Path('/Users/Alice/Documents/[Link]') # we can even use linux␣
↪like syntax on windows
print(todo_path.parts) # all parts of the path
print(todo_path.parent) # what directory is this file in?
print(todo_path.name) # what is the local name of the file?
print(todo_path.suffix) # what is the file extension?
print(todo_path.stem) # file without extension
# ... and many more

('C:\\', 'Users', 'Alice', 'Documents', '[Link]')


C:\Users\Alice\Documents
[Link]
.txt
todo

We can use these tools to create a path to another file in the same directory:

# create the new path by going to the parent and then adding a string␣
↪onto the path using the / (division) operator
# # this will combine the path with the string using the correct␣
↪separator for the platform
todo_path = Path(todo_path) # ignore this... it is necessary because it␣
↪was a pure path
other_file = todo_path.parent / "todo_backup.txt"
print(other_file)
# a path does not need to point to an existing file
print(other_file.exists()) # check if the file exists

5.2. Working with the File System 70


Python for Engineers

C:/Users/Alice/Documents/todo_backup.txt
False

Paths do not need to start at the root of the file system. They can also be relative to something like
the current working directory. The current working directory is usually where program is started
from or in this case it is the directory where VS Code is started from. We can get a path to the
current working directory by using the static [Link]() method.

[Link]() # get the current working directory

PosixPath('/absolute/path/to/working/directory')

The notes directory is in the current working directory, so we can also use a relative path to refer
to it.

relative_path = Path("notes")
print(relative_path.resolve()) # resolve the relative path to an␣
↪absolute path # BUILD COMMENT
# this is equivalent to
[Link]() / relative_path # we can also concatenate a path to another␣
↪path

PosixPath('/absolute/path/to/working/directory/notes')

If a path refers to a directory we can list all children of that directory with the iterdir method.

notes_dir = Path("notes")
print(notes_dir.is_dir()) # check if the path is a directory
print([Link]())
print(list(notes_dir.iterdir()))

True
/absolute/path/to/working/directory
[PosixPath('notes/math'), PosixPath('notes/mechanics'), PosixPath(
↪'notes/[Link]')]

To create the table of contents we will need to search for markdown files in the specified directory
and its subdirectories and sub subdirectories and so on. This is called recursively traversing the
file system which is also a good hint at how we can achieve this. Let’s write a function that takes a
path to a directory and lists its children, recursively calling itself if one of the children is a directory.

def recursive_traverse(directory: Path):


for child in [Link]():
# print(child) # print all files and directories
if [Link] == ".md": # or only print markdown files
print(child)
if child.is_dir(): # recurse if the path is a directory
recursive_traverse(child)

recursive_traverse(notes_dir)

notes/mechanics/mechanics 1/[Link]
notes/mechanics/mechanics 1/mechanics [Link]
(continues on next page)

5.2. Working with the File System 71


Python for Engineers

(continued from previous page)


notes/[Link]
notes/math/[Link]

Since this kind of thing is such a common pattern, there are already solutions for this in the stan-
dard library. However, it is still a useful concept to know for when these functions don’t serve your
specific needs.

# list directories and their contents


for directory, subdirs, files in notes_dir.walk():
print(f"{directory = }, {subdirs = }, {files = }")

directory = PosixPath('notes'), subdirs = ['mechanics', 'math'], files␣


↪= ['[Link]']
directory = PosixPath('notes/mechanics'), subdirs = ['mechanics 1'],␣
↪files = []
directory = PosixPath('notes/mechanics/mechanics 1'), subdirs = [],␣
↪files = ['[Link]', 'mechanics [Link]']
directory = PosixPath('notes/math'), subdirs = [], files = ['not␣
↪[Link]', '[Link]']

# search for all paths matching a "glob" pattern


# find all markdown files in the directory and nested directories
for file in notes_dir.rglob("*.md"): # does the same as notes_dir.glob(
↪"**/*.md")
print(file)

notes/[Link]
notes/math/[Link]
notes/mechanics/mechanics 1/[Link]
notes/mechanics/mechanics 1/mechanics [Link]

The Path class and the pathlib module have many more useful methods, if you need something
it probably exists, just look at the documentation.
Checkup Question
Consider the following code snippet:

file_path = Path("notes") / "[Link]"


output_path = file_path.parent / "[Link]"
print(output_path)

What will be the output?


• A: notes/[Link] (or notes\[Link] on Windows)
• B: notes/[Link]/[Link] (or notes\[Link]\[Link] on Win-
dows)
• C: [Link]
• D: An error because Path objects cannot be combined with strings using / for summary.
txt.
The correct answer is A: notes/[Link] (or notes\[Link] on Windows). The
.parent attribute file_path is a path to the notes directory. Using the / operator will then

5.2. Working with the File System 72


Python for Engineers

join the base path with [Link], resulting in notes/[Link].

5.3 Reading and Writing Files

To find the headings in the markdown files, we need to be able to look at their contents, meaning
we need to read them. We can do this by using the open function. It takes a path and gives us a
file object that we can read from. The object remembers our position in the file for the next read
operation.

# open the file in a with statement so that it is automatically closed


# at the end of the with block
markdown_file = Path("notes/math/[Link]")
with open(markdown_file) as f:
print([Link]()) # read all of the file
[Link](0) # go back to the beginning of the file
print([Link]()) # read a single line
print(list(f)) # read line by line
# [Link]("this does not work") # because the file is opened in read␣
↪mode

# My Math Notes
not part of a heading
#also not part of a heading

## Linear Algebra
### Matrices
## Differential Equations

# My Math Notes

['not part of a heading\n', '#also not part of a heading\n', '\n', '##␣


↪Linear Algebra\n', '### Matrices\n', '## Differential Equations\n']

By default, files are opened in read mode, if we want to write to a file we need to open it for writing
by specifying the mode as 'w' for write, by default this will create the file if it does not exist. There
are also other modes like 'a' for appending to a file or 'r+' for reading and writing.

table_of_contests_file = Path("[Link]")
print("file exists", table_of_contests_file.exists())

with open(table_of_contests_file, mode="w") as f:


# write a single line, don't forget the newline character
[Link]("# Table of Contents\n")
[Link](["- first entry\n", "- second entry\n"]) # write␣
↪multiple lines

print("file exists", table_of_contests_file.exists())


table_of_contests_file.unlink() # delete the file

file exists False


file exists True

The open function can treat a file as one of two different kinds: as a text file which it assumes by
default or as a binary file. Every file is technically a binary file, but files that only contain printable

5.3. Reading and Writing Files 73


Python for Engineers

characters can also be considered text files.


Checkup Question:
Which of these four snippets is the preferred way to read the lines of a file to a list?
• A:
file = open("[Link]", "r")
lines = [Link]()
[Link]()

• B:
with open("[Link]") as file:
lines = [Link]()

• C:
file = open("[Link]")
lines = []
for line in file:
[Link](line)

• D:
lines = open("[Link]", "r").read().splitlines()

The correct answer is B. The with statement ensures that the file is properly closed after its suite
finishes, even if an exception is raised. There is also no need to use a for loop to read the lines
as the method readlines already takes care of this.

5.4 Working with Strings

Finally, to be able to parse the contents of the files and generate the table of contents file, we still
need to learn a few things about working with strings. For example, we will need to find out if a
given line is level one or two heading.
For a line to be considered a heading it needs to start with a variable number of # characters
depending on the level of the heading, then at least one space and finally the actual heading text.
with open(markdown_file) as f:
lines = [Link]()

raw_headings: list[str] = []
for line in lines:
# strip removes leading and trailing whitespace, '\n' is whitespace
print(f'"{[Link]()}"', end="") # end='' prevents print from␣
↪adding a newline
if line[:2] == "# ": # are the first two characters '# '?
# easy to mess up, will always be false if the
# slice length and the string length are not the same
print(" -> top level heading")
raw_headings.append(line)
elif [Link]("## "): # use startswith instead of slicing,␣
(continues on next page)

5.4. Working with Strings 74


Python for Engineers

(continued from previous page)


↪harder to mess up
print(" -> second level heading")
raw_headings.append(line)
else:
print(" -> normal line")

"# My Math Notes" -> top level heading


"not part of a heading" -> normal line
"#also not part of a heading" -> normal line
"" -> normal line
"## Linear Algebra" -> second level heading
"### Matrices" -> normal line
"## Differential Equations" -> second level heading

Now that we identified what is and is not a heading, we still need to isolate the heading part. For
this we need to split the string on the first whitespace character. We also need to get in into the
correct form for a Markdown link, which means replacing the spaces with dashes and converting
the whole thing to lowercase. So given the raw heading string ## My Math Notes we want
somehow turn it into my-math-notes and determine that it is a level two heading.

for heading in raw_headings:


prefix, content = [Link](" ", 1) # split at the first space
content = [Link]() # strip potential leading white space and␣
↪the newline
heading_id = [Link](" ", "-").lower()
# markdown link [link text](link target)
# <> are needed when the link target contains spaces
# link a specific heading by using a '#' and the id of the heading to␣
↪link to
link = f"[{content}](<{file.as_posix()}#{heading_id}>)"
print(f'prefix: "{prefix}" content: "{content}"\nlink: "{link}"')

prefix: "#" content: "My Math Notes"


link: "[My Math Notes](<notes/mechanics/mechanics 1/mechanics [Link]#my-
↪math-notes>)"
prefix: "##" content: "Linear Algebra"
link: "[Linear Algebra](<notes/mechanics/mechanics 1/mechanics [Link]
↪#linear-algebra>)"
prefix: "##" content: "Differential Equations"
link: "[Differential Equations](<notes/mechanics/mechanics 1/mechanics␣
↪[Link]#differential-equations>)"

Using these basic string operations to check if a string fulfills specific criteria and then even more
operations to isolate certain parts can quickly become very complex and lead to hard to maintain
and most likely buggy code. If you encounter a situation like this, you should consider using
regular expressions. Regular expressions or regex are a special purpose language for exactly
this purpose and all major programming languages have support for them. You are not expected
to learn regex from this lecture, but you should know that it is a thing. To illustrate the basic
principles of regular expressions we will redo our example using regular expressions.
Regular expression can get very hard to read, so tools like regex101 and RegExr can be very
helpful in understanding, writing and testing them.

5.4. Working with Strings 75


Python for Engineers

import re

with open(Path("notes/math/[Link]")) as f:
whole_file = [Link]()

pattern = [Link](r"^(#{1,2})\s+(.+)$", [Link])


# `^` matches the beginning of a line
# `(...)` defines a capture group
# `#{1,2}` matches 1 or 2 '#' characters
# `\s+` matches one or more whitespace characters
# `(.+)` matches one or more of any character in another capture group
# `$` matches the end of a line
# `[Link]` means the input string is treated as multiple lines to␣
↪match individually

# iterate through all matches assigning the two capture groups to␣
↪`prefix` and `content`
for prefix, content in [Link](pattern, whole_file):
print(prefix, content)

# My Math Notes
## Linear Algebra
## Differential Equations

Play around with this regular expression and step through an interactive explanation on regex101.

5.5 String Interpolation with f-Strings

Another useful tool for working with strings that you should know about is string interpolation and
f-strings.
Whether we want to print out values for debug purposes or write them to a file in a structured way,
we often need to transform our numeric or other data into string form. While we could do all of this
by converting the individual to strings and then concatenating everything together with even more
strings, there is a much more elegant solution: string interpolation. There are multiple ways to
do this in Python, but we will be using f-strings. Like the name suggests, f-strings start with an f
before the quotes. Within the f-string, we can use curly braces to embed variables or expressions
into the string.

x = 2
f"i am a format string showing the number {x}"

'i am a format string showing the number 2'

print(f"multiple substitutions and even expressions: {x} and {x*x}")


print(f"use double braces for {{ no placeholder }}") # escape curly␣
↪braces with double braces

multiple substitutions and even expressions: 2 and 4


use double braces for { no placeholder }

5.5. String Interpolation with f-Strings 76


Python for Engineers

from datetime import datetime


import numpy as np

# f-strings also allow format specifiers that determine how the value is␣
↪converted to a string
print(f"floats: {[Link]}, with less decimals {[Link]:.4f}")
print(f"and in scientific notation: {[Link]:.2e}")
print(f"{42:04d}") # pad left with zeros
print(f"0x{42:x}, 0x{42:X} 0b{42:b}") # hexadecimals and binary
print(f"{x = }, {x**2 = }") # print an expression and its value (the␣
↪spaces matter)
print(f"current date and time: {[Link]():%Y-%m-%d %H:%M:%S}") #␣
↪format datetime objects
# and much more

floats: 3.141592653589793, with less decimals 3.1416


and in scientific notation: 3.14e+00
0042
0x2a, 0x2A 0b101010
x = 2, x**2 = 4
current date and time: 2026-04-13 20:34:58

Checkup Question
Consider the following code:

line = " # Chapter One: Introduction \n"


content = [Link]().split(":")[0].replace("# ", "")
print(f"'{content}'")

What will be the output?


• A: ' Chapter One'
• B: 'Chapter One'
• C: '# Chapter One'
• D: 'Chapter One: Introduction'
The correct answer is B: 'Chapter One'. The strip method removes leading and trailing
whitespace, the split method splits the string at the first colon, the [0] selects the first segment
of the split string, and the replace method removes the leading # from the string. The final
output is 'Chapter One' without any leading or trailing whitespace.

5.6 Putting it all together

Now that we figured out all the subproblems that we needed to solve our task, let’s put it all together
and write the program that generates the table of contents file. Starting with a docstring that
describes what the function should do.

from [Link] import Markdown, Pretty

def generate_table_of_contents(root_dir: Path, output_file: Path) ->␣


(continues on next page)

5.6. Putting it all together 77


Python for Engineers

(continued from previous page)


↪Markdown:
"""Recursively search for all markdown files in the root directory␣
↪and generate a table of contents

with all the headings and second level headings found in the files␣
↪and create hyperlinks to them.

Args:
root_dir (Path): The root directory to search for markdown files
output_file (Path): The file to write the table of contents to

Returns:
Markdown: The table of contents as a :class:`Markdown` object
"""

spaces = " "


extensions = {".md", ".MD", ".markdown"} # there are many more␣
↪common markdown extensions
result_lines = ["# Table of Contents", ""]

for directory, _, files in root_dir.walk():


# table of contents entry for directory
result_lines.append(spaces * len([Link]) + f" -
↪{[Link]}/")

for file_str in files:


file = directory / file_str
if [Link] not in extensions:
continue

# table of contents entry for file


result_lines.append(spaces * len([Link]) + f" - [
↪{[Link]}](<{file.as_posix()}>)")

# find and parse all headings


with open(file) as f:
pattern = [Link](r"^(#{1,2})\s+(.+)$", [Link])
headings: list[tuple[str, str]] = [Link](pattern, f.
↪read())

# generate headings with correct indentation


for level, heading in headings:
# table of contents entry for heading
link = f"{file.as_posix()}#{[Link](' ', '-').
↪lower()}"
result_lines.append(spaces * (len([Link]) +␣
↪len(level)) + f" - [{heading}](<{link}>)")

toc_contents = "\n".join(result_lines) + "\n"


with open(output_file, "w") as f:
[Link](toc_contents)

generate_table_of_contents(Path("notes"), Path("[Link]"))
Pretty(filename="[Link]")

5.6. Putting it all together 78


Python for Engineers

# Table of Contents

- notes/
- [[Link]](<notes/[Link]>)
- [Important Stuff](<notes/[Link]#important-stuff>)
- [How to Sign up for Exams](<notes/[Link]#how-to-sign-
↪up-for-exams>)
- mechanics/
- mechanics 1/
- [[Link]](<notes/mechanics/mechanics 1/[Link]>)
- [Exercises](<notes/mechanics/mechanics 1/[Link]
↪#exercises>)
- [mechanics [Link]](<notes/mechanics/mechanics 1/mechanics [Link]>)
- [Mechanics 1](<notes/mechanics/mechanics 1/mechanics [Link]
↪#mechanics-1>)
- [Statics](<notes/mechanics/mechanics 1/mechanics [Link]
↪#statics>)
- math/
- [[Link]](<notes/math/[Link]>)
- [My Math Notes](<notes/math/[Link]#my-math-notes>)
- [Linear Algebra](<notes/math/[Link]#linear-algebra>)
- [Differential Equations](<notes/math/[Link]#differential-
↪equations>)

To make sure that the Markdown is correct, we can also display the rendered Markdown:

# ... which will lead to problems when generating a pdf with latex
# Markdown(filename="[Link]") # rendered markdown

5.7 Summary

In this lecture we learned about:


• How to use pathlibs Path class to explore and manipulate the file system.
• How we can read and write the contents of individual files.
• Some common string manipulation tricks to create, manipulate and extract information from
strings.
• Using simple regular expressions to match and extract information from strings.
• How to combine all of the above to automate simple file system tasks.

5.7. Summary 79
CHAPTER

SIX

OBJECT-ORIENTED PROGRAMMING

Open in JupyterLite: Solved Blank


As a multi-paradigm language, Python also offers object-oriented programming features. But what
actually is object-oriented programming and why should we care? As the name implies, it evolves
around objects. These objects are basically containers for data that have functions attached to
them that act on that data and facilitate interaction between the outside world. We create these
objects from blueprints which we call classes that define what data fields every object of the
class should have and what operations called methods we can perform on it. When we create
an object from a certain class, we would call that object an instance of that class.
In previous lectures you already learned about built in classes like list and dict or custom
classes like numpy’s ndarray. Now it’s time to learn how to create classes of our own.

6.1 Defining a Class

We will do this by creating a class that should model complex numbers. Classes are defined using
the class keyword followed by the name of the class. Adhering to the PEP 8 style guide, we write
class names in PascalCase, e.g. capitalize every word and remove spaces.
class Complex:
# like functions, you should also put docstrings in classes to␣
↪document them

"""A class to represent complex numbers."""


# like type annotations in functions, these type annotations are also␣
↪optional
# and only serve the reader or your code editor to give hints
real_part: float
imaginary_part: float

# define methods inside (indented) the class


# member functions always (almost) take the instance to operate on as␣
↪the first argument
# the argument named 'self' by convention
def magnitude(self) -> float:
return (self.real_part**2 + self.imaginary_part**2) ** 0.5

z = Complex() # create an instance z, by 'calling' the class

# access attributes using dot notation


z.real_part = 3.0
z.imaginary_part = -4.0
(continues on next page)

80
Python for Engineers

(continued from previous page)

print(z) # only shows us what type of object and where it is


print(z.real_part, z.imaginary_part)
# explicitly get magnitude method from the class and pass the object to it
# print("magnitude:", [Link](z))
print("magnitude:", [Link]()) # implicitly call through the object,␣
↪the prefered way

<__main__.Complex object at 0x7b2655b044d0>


3.0 -4.0
magnitude: 5.0

In Python, objects are a bit like dictionaries, where the field names map to their values, only that
the syntax is a bit different.
In fact, they actually are dictionaries with some syntactic sugar on top. We can access the under-
lying dictionary using the special __dict__ attribute of the object.

z.__dict__

{'real_part': 3.0, 'imaginary_part': -4.0}

But where is the magnitude function that we defined? The magnitude function is not part of the
object itself, it belongs to the class so that’s where we find it.

Complex.__dict__

mappingproxy({'__module__': '__main__',
'__annotations__': {'real_part': float, 'imaginary_part':␣
↪float},
'__doc__': 'A class to represent complex numbers.',
'magnitude': <function __main__.[Link](self) ->
↪ float>,
'__dict__': <attribute '__dict__' of 'Complex' objects>,
'__weakref__': <attribute '__weakref__' of 'Complex'␣
↪objects>})

6.2 Special Methods

Creating objects like this is a bit cumbersome and limiting. To make this a bit more ergonomic we
can define the special __init__ method. This method is called with the parameters we pass to
the “class call” and is used to initialize the object. This is often also called a constructor, because
it constructs the object.

class Complex:
real_part: float
imaginary_part: float

def __init__(self, real: float, imag: float):


self.real_part = real
self.imaginary_part = imag
(continues on next page)

6.2. Special Methods 81


Python for Engineers

(continued from previous page)

# another special method to make the class print out more usefull
# the method should be defined, such that `eval(repr(obj)) == obj`
# eval takes a string and executes it as python code returning the␣
↪result
def __repr__(self) -> str: # short for representation
return f"Complex({self.real_part}, {self.imaginary_part})"

# python will call the __add__ method when we write Complex(a, b) +␣


↪Complex(c, d)
# this is called operator overloading and is available for most␣
↪operators
def __add__(self, other: "Complex") -> "Complex":
return Complex(self.real_part + other.real_part, self.imaginary_
↪part + other.imaginary_part)

z = Complex(1.0, 2.0)
print(z)
# print([Link](Complex(3.0, 4.0)))
print(z + Complex(3.0, 4.0))
print(f"{repr(z)==repr(eval(repr(z))) = }") # compare the string since␣
↪we didn't define equality

Complex(1.0, 2.0)
Complex(4.0, 6.0)
repr(z)==repr(eval(repr(z))) = True

These methods are also called dunder methods because they are surrounded by double un-
derscores and there are many more of them. As we already saw with __dict__ there are also
dunder attributes. Unless you either know what you are doing or want to implement one of
these special methods you should not start you own method or attributes names with a double
underscore.

6.3 Complex Numbers in Python

This was only made to serve as an example, you should not reimplement complex numbers your-
self. Complex numbers are actually built into Python and you can create them by using the imag-
inary unit j or their constructor complex.

w = 2 + 3j
print(w)
print([Link], [Link])
print(type(w))
print(w + 2)

(2+3j)
2.0 3.0
<class 'complex'>
(4+3j)

6.3. Complex Numbers in Python 82


Python for Engineers

6.4 Encapsulation

Another important concept of object-oriented programming is encapsulation. As classes become


more complicated, it gets easier to make mistakes and put them in an inconsistent state. To
prevent this, we might want to encapsulate the data and only allow changes to the object through
methods where we can make sure that the object is always left in a consistent state.
For example, coordinates on a sphere should always have a latitude between -90 and 90 degrees
and a longitude between -180 and 180 degrees. So let’s implement a class that ensures this.

from math import radians

class SphereCoordinate:
# the leading underscore tell people that this is a private attribute
# and should not be accessed from outside the class
_lat_deg: float
_lon_deg: float

def __init__(self, lat: float, lon: float):


if not (-90 <= lat <= 90):
raise ValueError(f"Latitude {lat} out of range, must between -
↪90 and 90")
if not (-180 <= lon <= 180):
raise ValueError(f"Longitude {lon} out of range, must between␣
↪-180 and 180")
self._lat_deg = lat
self._lon_deg = lon

def __repr__(self) -> str:


return f"Coordinate({self._lat_deg}, {self._lon_deg})"

# this is a decorator, for now you can think of them as special␣


↪anotations
# that tell python to do something special with the thing below them
# you should only use property for very simple and cheap to execute␣
↪functions

# this allows us to write coord.lat_deg instead of coord.lat_deg()


@property
def lat_deg(self) -> float:
return self._lat_deg

# this is called when assigning to the attribute eg. coord.lat_deg =␣


↪10
@lat_deg.setter
def lat_deg(self, lat: float):
if not (-90 <= lat <= 90):
raise ValueError(f"Latitude {lat} out of range, must between -
↪90 and 90")

self._lat_deg = lat

@property
def lon_deg(self) -> float:
return self._lon_deg

# they also allow us to make attribute like things that are computed␣
↪from other
(continues on next page)

6.4. Encapsulation 83
Python for Engineers

(continued from previous page)


# attributes, so they don't need to be stored and therefore cant be␣
↪inconsistent
@property
def lat_rad(self) -> float:
return radians(self._lat_deg)

# ...

coord = SphereCoordinate(52.38, 9.72)


# bad_coord = SphereCoordinate(100, 200) # this will raise an error
print(coord)

coord.lat_deg = 50 # ok, valid value


print(coord.lat_deg) # use the getter method
# coord.lat_deg = 100 # this will raise an exception because the value␣
↪is out of range
# coord._lat_deg = 100 # this works but the underscore tells us that we␣
↪shouldn't do this

print(coord.lat_rad) # automatically computes the angle in radians

Coordinate(52.38, 9.72)
50
0.8726646259971648

Encapsulation does not only help in ensuring consistent state, it also isolates the interface from the
implementation. For example, if for some reason we want to change the internal representation
of the coordinates from degrees to radians, this implementation would allow us to do this since all
access goes through the getters and setter where we can make the necessary conversions and
consistency checks.
One thing to keep in mind is that there is nothing to be gained by implementing a getter and setter
which read or write to the underlying field without any additional logic. It adds nothing and behaves
just like a public field. Properties are most useful for read only or computed fields. Your first line
of defense against inconsistency should not be a setter which validates the input, it should be to
disallow any change at all since it is often unnecessary and fewer changing things generally lead
to less confusion.

6.5 Inheritance

Like other languages which offer OOP features, Python, of course, also supports inheritance. We
will only briefly go over this, as you are not going to use it as much as people who (religiously)
teach OOP would have you believe.
Inheritance is a way to model relationships between types and lets them reuse common function-
ality. Inheritance models an is a relationship, like an engineering student is a student or a student
is a person. In this case, the engineering student would inherit from the student and the student
from the person. By deriving / inheriting from another class, the subclass inherits all of the fields
and methods of the parent class and can add additional fields and methods or override existing
ones. In Python, this would look like this:

6.5. Inheritance 84
Python for Engineers

class Person: pass # use pass to define an empty class


class Student(Person): pass # inherit all attributes and methods from␣
↪parent class in parentheses
class Professor(Person): pass
class EngineeringStudent(Student): pass

# to check wether a *is a* B we use isinstance(a, B)


# isinstance(a, B) is similar to type(a) == B but it also checks if a is␣
↪instance of a subclass of B
print("person is a Person", isinstance(Person(), Person)) # Person is a␣
↪Person
print("student is a Student", isinstance(Student(), Person)) # Student␣
↪is a Person
# type(student) == Student is false since Student is a subclass of Person
print("type(student) is Person", type(Student()) is Person)
print("person is a Student", isinstance(Person(), Student)) # generic␣
↪Person is not a Student
# *is a* is transitive
print("engineering_student is a Person", isinstance(EngineeringStudent(),␣
↪Person))

person is a Person True


student is a Student True
type(student) is Person False
person is a Student False
engineering_student is a Person True

Since is a is transitive, an EngineeringStudent is of course also a Person. This means


that everything operation or function that takes a Person should also work with an Engineer-
ingStudent without any changes. This is called the Liskov substitution principle.
Using inheritance, we can put our data and functions that all people have in common into the
Person class, and then put special additional functionality for a student, like for example the
grades, in the Student class. By inheriting from Person, Students have all the functionality and
data from a Person like name and age.
Let’s program part of that example:

class Person:
def __init__(self, name: str, age: int):
[Link] = name
[Link] = age

def say_hello(self):
print(f"Hello, my name is {[Link]} and I am {[Link]} years␣
↪old.")

def __repr__(self) -> str:


return f"Person({[Link], [Link]})"

class Student(Person):
def __init__(self, name: str, age: int, grades: dict):
Person.__init__(self, name, age) # delegate to the parent classes
[Link] = grades

# since we inherit from Person we get say_hello for free


(continues on next page)

6.5. Inheritance 85
Python for Engineers

(continued from previous page)

def __repr__(self) -> str: # override the __repr__ function with one␣
↪that includes grades
return f"Student({[Link], [Link], [Link]})"

people = [Person("Alice", 22), Student("Bob", 24, {"Maths": 1.7,


↪"Mechanics": 2.0})]

for person in people:


print("type:", type(person))
print(person)
person.say_hello()
print("")

bob = people[1]
print(bob)
print("is Student?", isinstance(bob, Student))
print("is Person?", isinstance(bob, Person))

type: <class '__main__.Person'>


Person(('Alice', 22))
Hello, my name is Alice and I am 22 years old.

type: <class '__main__.Student'>


Student(('Bob', 24, {'Maths': 1.7, 'Mechanics': 2.0}))
Hello, my name is Bob and I am 24 years old.

Student(('Bob', 24, {'Maths': 1.7, 'Mechanics': 2.0}))


is Student? True
is Person? True

This is a typical “academic” example of inheritance and should be taken as such. In practice you
would likely not implement such a thing like this. Especially in multi-paradigm languages where
there are many other tools that can solve the sample problems. You rarely actually need to use
inheritance and most of the time you are better off just putting the parent class into the would-be
child class, which is called composition. See here for the reasoning behind this.
Real world example of inheritance are often much more abstract and don’t have real world analo-
gies like the one above. One example that you already learned about is [Link] and
its related classes like PosixPath, WindowsPath and PurePath. This has little to do with the
common real world comparison of literal inheritance. Here we can also see the Liskov substitution
principle in action: In lecture 4 we defined all our function to operate on Path objects even though
they actually received PosixPath or WindowsPath objects.

6.5. Inheritance 86
Python for Engineers

6.6 Dataclasses & Namedtuples

In a lot of cases, you just want to have a simple data container without any complicated initialization
functions or other validation. In this case you should consider using a dataclass where the
__init__ and __repr__ functions are automatically generated for you. This make the code
much simpler and is one less opportunity to make a mistake.

from dataclasses import dataclass

@dataclass # this decorator it will modify the class


class Point2D:
# no init function, we just declare the members at class level
x: int # no default value for x and y
y: int
# if we provide a default value it will become an optional argument␣
↪in the constructor
name: str = "default"

# you can still define your own methods


def norm(self) -> float:
return (self.x**2 + self.y**2) ** 0.5

p = Point2D(1, 2) # the init function is created automatically


print(p) # the repr function is also created automatically
p_named = Point2D(3, 4, "P") # give a name
print(p_named, p_named.norm())

Point2D(x=1, y=2, name='default')


Point2D(x=3, y=4, name='P') 5.0

But to the minimalist Python programmer, even this might be too much, and they would really just
want a fancier tuple with named fields. This is where namedtuple comes in, which is a function
that generates simple tuple like classes.

from collections import namedtuple

Person = namedtuple("Person", ["first_name", "last_name", "birth_year"])

guido = Person("Guido", "van Rossum", 1956)


print(guido)
print(guido.first_name)

# named tuples are actual tuple so they inherit all of tuple's properties
print(guido[0]) # like indexing
first, last, year = guido # and destructuring
print(first, last, year)
# and being immutable
# guido.first_name = "Bob" # this will raise an error

Person(first_name='Guido', last_name='van Rossum', birth_year=1956)


Guido
Guido
Guido van Rossum 1956

6.6. Dataclasses & Namedtuples 87


Python for Engineers

6.7 Summary

In This lecture we learned about:


• The difference between objects and classes and how to define them.
• Special methods like __init__ and __repr__ and what the interpreter does with them.
• Object-oriented concepts like encapsulation and inheritance.
• Utilities like dataclass and namedtuple for simple composite types.

6.8 Appendix - Decorators

Unlike their usage might suggest, decorators aren’t magic. To get a better understanding of how
they work, we will implement one ourselves.
In essence, a decorator is just a function that takes another function and wraps or modifies it in
some way, or in the case of a class it takes the whole class and modifies it.
Let’s implement a decorator that prints out debug information when calling the decorated function:

from typing import Callable, Any


from functools import wraps

# the decorator takes a callable taking one argument and returns a␣


↪callable

def unary_debug_print(func: Callable[[Any], Any]) -> Callable[[Any], Any]:


# create a wrapper function that forwards the parameter and result␣
↪and print out debug info
# @wraps(func) # preserve the name and docstring of the original␣
↪function
def wrapper(arg: Any) -> Any:
print(f"Calling {func.__name__} with argument: {arg}")
result = func(arg)
print(f"{func.__name__} returned: {result}")
return result

# wrapper.__name__ = func.__name__ # we can patch the name and other␣


↪dander attributes
# to hide the fact that the wrapper is not the original function, or␣
↪just use @wraps decorator

return wrapper # the function we decorate will be replaced by the␣


↪wrapper

To apply a function as a decorator we use the @ syntax before the function’s definition. The dec-
orated function is replaced by whatever the decorator function returns:

@unary_debug_print
def square(x: float) -> float:
return x**2

res = square(2.0) # we get the debug output implemented in the decorator

(continues on next page)

6.7. Summary 88
Python for Engineers

(continued from previous page)


print(square.__name__) # prints `wrapper` because the `square` function␣
↪was replaced by the wrapper

Calling square with argument: 2.0


square returned: 4.0
wrapper

The @ decorator syntax is just syntactic sugar to make decorators more convenient to use. We
can achieve the same results by “decorating” the function manually:

def cube(x: float) -> float:


return x**3

cube = unary_debug_print(cube) # this does the same as the add syntax

res = cube(3.0) # we get the debug output implemented in the decorator

Calling cube with argument: 3.0


cube returned: 27.0

6.8. Appendix - Decorators 89


CHAPTER

SEVEN

BITS AND PIECES

Open in JupyterLite: Solved Blank


In this lecture, we are going to talk about some Python features that, while not strictly necessary to
write functioning code, can make your life a lot easier. You don’t have to remember all the details
now, but consider revisiting this from time to time to improve your coding style.

7.1 Unpacking Assignments

We can assign the elements of a tuple, list, or any other iterable to individual variables using
unpacking assignments.

# simple unpacking
point = (1, 2) # bracket are optional
x, y = point
print(x, y)

# we can also do this with nested data structures by matching the␣


↪structure on the assignment side
named_point = ("p1", (3, 4))
name, (x, y) = named_point
print(name, x, y)

# since we can create the duple on the right side, we can also use this␣
↪to do multiple assignments at one
# this can make for nice looking code
a, b = 5, 6
print(a, b)
# we can also use this to swap variables
a, b = b, a
print(a, b)

1 2
p1 3 4
5 6
6 5

This can also be used to intuitively return multiple values from a function:

def sum_and_product(a, b):


sum = a + b
(continues on next page)

90
Python for Engineers

(continued from previous page)


prod = a * b
return sum, prod

s, p = sum_and_product(2, 3)
print(s, p)

5 6

7.2 Better Iteration with zip() and enumerate()

There are many situations where we want to iterate over both the indices and the values of a
list. Or we want to iterate over two or more iterables at the same time. The naive solution to this
problem is to just iterate over the indices and then use square brackets to access the elements.
This is not very elegant and is prone to errors. Instead, enumerate() and zip() should be
used to write more concise code.

names = ["Alice", "Bob", "Charlie"]


ages = [22, 19, 23]

# DON'T DO THIS
for i in range(len(names)):
print(f"{names[i]} is {ages[i]} years old")
print("")

# do this instead
for name, age in zip(names, ages, strict=True):
# use strict=True so that the zip will throw an error if the lists␣
↪are not the same length
print(f"{name} is {age} years old")

Alice is 22 years old


Bob is 19 years old
Charlie is 23 years old

Alice is 22 years old


Bob is 19 years old
Charlie is 23 years old

In other situations, you might also need the index of the items while iterating over them. Instead
of iterating over just the index and then using the [] operator to access the element, you should
use enumerate().

# DON'T DO THIS
for i in range(len(names)):
print(f"{names[i]} has index {i}")
print("")

# do this instead
for i, name in enumerate(names):
print(f"{name} has index {i}")

7.2. Better Iteration with zip() and enumerate() 91


Python for Engineers

Alice has index 0


Bob has index 1
Charlie has index 2

Alice has index 0


Bob has index 1
Charlie has index 2

# we can also use these two together


for i, (name, age) in enumerate(zip(names, ages, strict=True)):
# nested unpacking assignment
print(f"{name} is {age} years old and has index {i}")
print("")

# zip can also take any number of arguments


semesters = [5, 1, 6]
for name, age, semester in zip(names, ages, semesters, strict=True):
print(f"{name} is {age} years old and studying in semester {semester}
↪")

Alice is 22 years old and has index 0


Bob is 19 years old and has index 1
Charlie is 23 years old and has index 2

Alice is 22 years old and studying in semester 5


Bob is 19 years old and studying in semester 1
Charlie is 23 years old and studying in semester 6

7.3 List Comprehensions

List comprehensions are a very powerful and elegant language feature in Python. They can be
used to transform and filter lists. The source can also be a generator expression or any other
iterable. With all these features combined, complex series can be generated with little code.

numbers = [0, 0.5, 1, 1.5, 2]


squared_numbers = [num**2 for num in numbers] # transform list
print(squared_numbers)
# the source can also be a generator (or other iterable)
squared_integers = [num**2 for num in range(10)]
print(squared_integers)

# filter list (or generator)


odd_squares = [num for num in squared_integers if num % 2 == 1]
print(odd_squares)

# mapping / transformation and filtering combined


even_powers_of_two = [2**i for i in range(10) if i % 2 == 0]
print(even_powers_of_two)

[0, 0.25, 1, 2.25, 4]


[0, 1, 4, 9, 16, 25, 36, 49, 64, 81]
[1, 9, 25, 49, 81]
[1, 4, 16, 64, 256]

7.3. List Comprehensions 92


Python for Engineers

List comprehensions can also be nested within each other to generate or transform multidimen-
sional lists. Just like regular for loops, they can also be combined with structured bindings to
iterate through multiple lists at once.

products = [[i * j for j in range(10)] for i in range(5)]


print(products)

[[0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 1, 2, 3, 4, 5, 6, 7, 8, 9], [0, 2,␣


↪4, 6, 8, 10, 12, 14, 16, 18], [0, 3, 6, 9, 12, 15, 18, 21, 24, 27],␣
↪[0, 4, 8, 12, 16, 20, 24, 28, 32, 36]]

The equivalent of list comprehensions also exists for dictionaries and sets, as dictionary and set
comprehensions.

my_set = {val for val in range(3)}


print(my_set)
my_dict = {val: str(val) for val in range(3)}
print(my_dict)

{0, 1, 2}
{0: '0', 1: '1', 2: '2'}

Tuple comprehensions do not exist, but something similar can be achieved by passing a generator
expression to the tuple() constructor.

# parentheses create a generator expression, not a tuple


generator = (num for num in range(3))
print(generator)

# but we can pass this generator to the tuple constructor, to get␣


↪something close to a tuple comprehension
my_tuple = tuple(num for num in range(3))
print(my_tuple)

<generator object <genexpr> at 0x74c95c23ad40>


(0, 1, 2)

7.4 Reference Semantics

In Python, assigning a variable to another name does not create a copy of the object, but only
a shared reference to it. This can be very confusing when working with mutable objects, so it is
important to know what’s really going on.
We can check if two variables refer to the same object by comparing their id()s or using the is
operator which does the same thing under the hood.

my_list = [1, 2, 3]
the_same_list = my_list
the_same_list[1] = 5
print(f"{my_list = }, {the_same_list = }") # both lists are changed
print(f"{id(my_list) == id(the_same_list) = }") # because they are the␣
↪same object

(continues on next page)

7.4. Reference Semantics 93


Python for Engineers

(continued from previous page)


print(f"{my_list is the_same_list = }") # you should use `is` for␣
↪identity comparison

my_list = [1, 5, 3], the_same_list = [1, 5, 3]


id(my_list) == id(the_same_list) = True
my_list is the_same_list = True

my_list = [1, 2, 3]
a_different_list = my_list.copy()
print(f"{my_list == a_different_list = }") # the values are the same␣
↪(equal)
print(f"{my_list is a_different_list = }") # but the objects are␣
↪different (not identical)
# so changing one doesn't change the other
a_different_list[1] = 5
print(f"{my_list = }, {a_different_list = }")

my_list == a_different_list = True


my_list is a_different_list = False
my_list = [1, 2, 3], a_different_list = [1, 5, 3]

These reference semantics also apply to function arguments. When we are passing a variable to
a function, we are only giving it a shared reference to the same data and the function works on
the original object.

my_list = [1, 2, 3]

def do_something(data: list) -> list:


data[1] = 42
return data

changed_list = do_something(my_list)
print(changed_list) # the function returns the changed list
print(my_list) # the original list is also changed
# this is because changed_list and my_list are both a reference to the␣
↪same object
print(my_list is changed_list)

[1, 42, 3]
[1, 42, 3]
True

Often times, this is fine, especially when performing read-only operations on the list or when it
doesn’t matter what happens to the list. As a general rule, it is not a good idea to modify input
arguments, unless it is obvious from the function name. If you need something like the input to
modify, just make a copy of it.

my_list = [1, 2, 3]

def do_something_else(data: list) -> list:


new_data = [Link]()
new_data[1] = 42
return new_data
(continues on next page)

7.4. Reference Semantics 94


Python for Engineers

(continued from previous page)

a_new_list = do_something_else(my_list)
print(my_list)
print(a_new_list)

[1, 2, 3]
[1, 42, 3]

This can be especially treacherous when working with mutable default arguments that are sup-
posed to be mutated in the function. The default argument is only evaluated once when the func-
tion is defined, and not every time the function is called as one might expect.

def append_one(data: list = []) -> list:


[Link](1)
return data

print(append_one([1, 2, 3])) # everything is fine when you provide a list


print(append_one()) # this will use the default list
print(append_one()) # subsequent calls will use the same (identical) list

[1, 2, 3, 1]
[1]
[1, 1]

To implement the desired behavior correctly, we can use None as the default value and check for
it in the function body.

def append_one_fixed(data: list | None = None) -> list:


data = data or [] # equivalent to `data if bool(data) else []`
[Link](1)
return data

print(append_one([1, 2, 3]))
print(append_one_fixed()) # new list is created in the function body
print(append_one_fixed())

[1, 2, 3, 1]
[1]
[1]

Another common pitfall is creating a nested list with the multiplication operator.

one_dim = [0] * 3
one_dim[1] = 1
print(one_dim) # works as expected

two_dim = [[0] * 3] * 3
two_dim[1][1] = 1
print(two_dim) # unintuitive behavior
print(two_dim[0] is two_dim[1]) # all the inner lists are the same␣
↪object (identical)

# use a list comprehension to create a new list for each row


two_dim_fixed = [[0] * 3 for _ in range(3)]
(continues on next page)

7.4. Reference Semantics 95


Python for Engineers

(continued from previous page)


two_dim_fixed[1][1] = 1
print(two_dim_fixed) # now everything works as expected

[0, 1, 0]
[[0, 1, 0], [0, 1, 0], [0, 1, 0]]
True
[[0, 0, 0], [0, 1, 0], [0, 0, 0]]

7.5 Lambdas

Sometimes we need to pass functions to other functions. These are often very simple, and using
the def syntax makes for unpleasant code. The solution to this is to use lambda expressions.

from typing import Callable, Any

# the following two are more or less equivalent


# my_func = lambda [params] : [expression]
# def my_func([params]):
# return [expression]

def accumulate(arr: list, func: Callable):


res = arr[0]
for val in arr[1:]:
res = func(val, res)
return res

numbers = [1, 4, 2, 3]
print("sum", accumulate(numbers, lambda x, y: x + y)) # lambda [params]␣
↪: [expression]
print("prod", accumulate(numbers, lambda x, y: x * y))
print("max", accumulate(numbers, max))
magic_number = 3
# lambdas can access the variables of the scope they are created in
print("bad idea", accumulate(numbers, lambda x, y: magic_number))

sum 10
prod 24
max 4
bad idea 3

Passing functions to other functions is part of the functional programming paradigm. This can
make for very elegant code, but in some situations, it is a lot harder to understand than the imper-
ative style programming.

7.5. Lambdas 96
Python for Engineers

7.6 Function Parameters

There are a lot of ways to specify and pass function parameters in Python. We can set default
values, pass parameters both by position or name, and we can use * and ** to pass multiple
parameters at once.

def f1(a, b, c): # the usual way


print(f"f1: {a=}, {b=}, {c=}")

f1(1, 2, 3) # the "normal" way


params = (4, 5, 6)
f1(*params) # unpack a tuple or list of params
kwargs = {"c": 7, "b": 8, "a": 9} # the order in the dict does no matter
f1(**kwargs) # unpack a dict of params
print("")

def f2(a=1, b=2, c=3): # with default values


print(f"f2: {a=}, {b=}, {c=}")

f2()
f2(4) # pass the first of the defaults
f2(b=9) # only specify a value for b

f1: a=1, b=2, c=3


f1: a=4, b=5, c=6
f1: a=9, b=8, c=7

f2: a=1, b=2, c=3


f2: a=4, b=2, c=3
f2: a=1, b=9, c=3

def f3(*args): # absorb all positional parameters into tuple args


print(f"f3: {args=}")

f3()
f3(1, 2, 3)
print("")

# the names "args" and "kwargs" are just conventions, the stars are what␣
↪matters
def f4(a, *args, **kwargs): # absorb all other keyword parameter into␣
↪dict kwargs

print(f"f4: {a=}, {args=}, {kwargs=}")

f4(1)
f4(1, 2, 3, d=4, e=5)

f3: args=()
f3: args=(1, 2, 3)

f4: a=1, args=(), kwargs={}


f4: a=1, args=(2, 3), kwargs={'d': 4, 'e': 5}

7.6. Function Parameters 97


Python for Engineers

# check if a specific keyword parameter was passed


def f5(**kwargs):
if "a" in kwargs:
print(f"a was passed with value {kwargs['a']}")
else:
print("a was not passed")

f5(a=2, b=3)
f5(c=4)

a was passed with value 2


a was not passed

We can also force the user of our function to pass parameters by position or value by including the
special / and * placeholder arguments in the function’s signature. All parameters before the /
must be passed by position, and all parameters after the * must be passed by name. Everything
in between can be passed by both.

def f6(pos_only, /, pos_or_kw, *, kw_only):


pass

f6(1, 2, kw_only=3) # ok
f6(1, pos_or_kw=2, kw_only=3) # also ok

f6(1, 2, 3) # error because `kw_only` must be passed as a keyword␣


↪argument

------------------------------------------------------------------------
↪---
TypeError Traceback (most recent call␣
↪last)
Cell In[22], line 1
----> 1 f6(1, 2, 3) # error because `kw_only` must be passed as a␣
↪keyword argument

TypeError: f6() takes 2 positional arguments but 3 were given

f6(pos_only=1, pos_or_kw=2, kw_only=3) # error because `pos_or_kw` must␣


↪be passed as a positional argument

------------------------------------------------------------------------
↪---
TypeError Traceback (most recent call␣
↪last)
Cell In[23], line 1
----> 1 f6(pos_only=1, pos_or_kw=2, kw_only=3) # error because `pos_or_
↪kw` must be passed as a positional argument

TypeError: f6() got some positional-only arguments passed as keyword␣


↪arguments: 'pos_only'

This is useful for forcing optional parameters to be passed only as keyword arguments, which
allows adding new parameters changing their order without potentially breaking the code that
relied on the order of the parameters.

7.6. Function Parameters 98


Python for Engineers

7.7 Code Style

Writing good code is not just about writing code that works, it is also about writing code that is
easy to read and understand. There is no one set of rules that you can follow to write better code,
but there are a few rules that, when applied, will generally lead to nicer code.
You will not get to apply all of these conventions, as the exercises are not complex enough for
some of them to be relevant. However, you should still be aware of them and try to apply them
when you can.

7.7.1 Naming Conventions

Phil Karlton once famously said: “There are two hard things in computer science: cache inval-
idation and naming things.” Whether true or not, naming things is certainly very important, and
sometimes it is really hard to find a good name for something. Nevertheless, there are a few rules
that you should follow when naming things in Python.
The first isn’t actually the names, but how you write them. They should follow a consistent style.
When this style is not already specified by the project you are working on, you should follow the
PEP 8 style guide. The essence of that is:
• Use snake_case for variable and function names
• Use PascalCase for class names
• Use UPPER_SNAKE_CASE for constants
• Don’t use single character variable names, except for counters in loops (we sometimes
break this rule when doing math)
• Don’t use cryptic abbreviations without vowels. It’s fine if your variable names consist of
multiple words, but don’t make them too long either. If you have a very big expression, you
can still use local shorthands to reduce visual clutter.
Try to name your variables and functions so that it is obvious what the variable contains / how it
is calculated, or what the function does without having to look up the definition. While comments
are a useful feature to give the reader information about things that are not immediately obvious
from threading the code, it is event better to have the code be self-explanatory. Instead of writing:

ll = 11 # length of the left list


lr = 42 # length of the right list

Use more descriptive variable names to make the comments unnecessary:

length_left_list = 11
length_right_list = 42

7.7. Code Style 99


Python for Engineers

7.7.2 Automatic Formatting & Linting

Luckily you do not have to carry all of this burden yourself. There are tools out there which can
help you keep you code tidy and warn you of common malpractices. Autoformatters like Black
for Python mostly add and remove whitespace to make everything look consistent. Linters like
Ruff on the other hand will warn you of common mistakes or anti-patterns and sometimes offer
automatic fixes.
After installing Black you will need to set these settings in VS-Code’s user or workspace
[Link] to configure it:

{
"[Link]": true,
"[Link]": true,
"[python]": {
"[Link]": "[Link]-formatter"
},
"[Link]": ["-l", "99"]
}

The setting [Link] tells VS-Code to format the code every time you save
the file. For this to work in notebooks, you also need to explicitly enable it with notebook.
formatOnSave. The last setting [Link] increases the allowed maximum
line length to 99 characters to make better use of modern monitors.
It is also possible to set these settings in the settings UI by just searching for the names of the
setting.

7.7.3 Use More Functions

Functions are probably the most powerful abstraction there is, use them! Once you find yourself
writing the same code twice, you should think about writing a function for it. If you write it three
times, you should have already written the function. Often, it will also make sense to write a
function even if you only use it once. You can do this to make it clear that this is a meaningful sub-
operation of your function. This also helps to prevent your functions from turning into a 100-line
messes.

7.7.4 Use Early Returns

You will often encounter situations where you need to check that multiple conditions are met before
doing the actual work. A bad way to do this is to keep nesting if statements.

def make_even(x: int | None) -> int | None:


if x is not None:
if x % 2 == 1:
# here the normal case is two levels deep
return 2 * x # image this line is a very complicated block␣
↪of code
else:
return x
return None

7.7. Code Style 100


Python for Engineers

As you can see, this will quickly lead to a lot of indentation and make the code hard to read.
Instead, you should use early returns. This means that you check for error cases first and return
early if you find one. This will make the code much easier to read.

def better_make_even(x: int | None) -> int | None:


if x is None:
return None
if x % 2 == 0:
return x

# here the normal case is at base indentation level


return 2 * x # image this line is a very complicated block of code

7.7.5 Use Type Hints

Use type hints! They make it easier to understand what a function does and how it is used. They
will also help your IDE give you better suggestions for auto-completion and find errors in your code.

7.7.6 Prioritize Readability

Whenever you are unsure about how to write something, you should (almost) always go for the
more readable option. This might not always be the most efficient one or the one with the fewest
lines of code, but it will likely save you and everyone else time. Complexity is a fairly good measure
of how hard it is to understand something, so it is important to always strive to write the least
complex code possible.

7.7.7 Don’t Be Afraid to Start Over

It is unlikely that you will get “it” right the first time. Burning it all down and starting from scratch
is often way less work than you think, since most of the work was probably trying to figure out the
problem, not writing the code. You can carry over that experience to the new implementation and
write it much faster and better without the holes that you dug yourself into on your first attempt.

7.8 Summary

In this lecture we learned about:


• Using unpacking assignments for more concise code and returning multiple values from
functions.
• Using iteration tools like zip() and enumerate() to write cleaner and less error-prone
loops.
• How to write list (and dict/set) comprehensions to generate or transform lists concisely.
• The reference semantics of Python and how to avoid common pitfalls.
• Lambda functions and what we can do with them.
• More advanced ways to pass and receive function parameters.

7.8. Summary 101


Python for Engineers

• Some guidelines on how to write better code.

7.8. Summary 102


CHAPTER

EIGHT

PANDAS

Open in JupyterLite: Solved Blank


Although Pandas is mostly known for its uses in data science and related fields, its two core data
structures the Series and the DataFrame are also very useful for general data manipulation
tasks. As with most of these established libraries, the documentation contains a lot of of useful
resources to get started.

import numpy as np
import pandas as pd
import [Link] as plt

[Link]("[Link]")

8.1 Series for Ordered Data

Up to this points we have already spent some time getting to know the NumPy ndarray and it
has served us well for most numerical needs. But just because we were able to manage, does
not mean we should not look for better tools.
Often times we have some NumPy array that contains our index points like times ts and some
data vals that represent the values at these time points. These variables belong together but
there really isn’t anything actually keeping them together and forcing them to be consistent. If they
only exist as short-lived temporary values this is mostly unproblematic, but if we use them more
and often pass them between functions they really should be contained in some data structure to
ensure consistency and make the code simpler and less error-prone.
This is where Panda’s Series comes in, it can hold the data and the time index in one object and
also implements a lot of useful methods for working with the data.

ts = [Link](0, 3600, 7) # time in seconds


temp_deg = [Link]([17, 18, 20, 21, 21, 19, 20]) # temperature in␣
↪degrees Celsius

# bundle the free variable time and the dependent variable temperature␣
↪into a pandas Series

temperature_series = [Link](temp_deg, index=ts)


[Link](temperature_series); # uses the index as the x-axis

103
Python for Engineers

# mathematical operations act on the series' values and return a new␣


↪series with the same index
# in many situations series behave just like numpy arrays
temp_kelvin = temperature_series + 273.15
temp_kelvin

0.0 290.15
600.0 291.15
1200.0 293.15
1800.0 294.15
2400.0 294.15
3000.0 292.15
3600.0 293.15
dtype: float64

8.2 DataFrames as Collections of Series

If we have another series of values that is defined for the same index, we can use a DataFrame,
it is essentially a dictionary of Series objects that share the same index.

temperatures = [Link](
{"inside_temp": temp_deg, "outside_temp": temp_deg[::-1] - 10},␣
↪index=[Link](ts, name="time")
)
temperatures # dataframes have nice html representations in notebooks

inside_temp outside_temp
time
0.0 17 10
600.0 18 9
1200.0 20 11
1800.0 21 11
(continues on next page)

8.2. DataFrames as Collections of Series 104


Python for Engineers

(continued from previous page)


2400.0 21 10
3000.0 19 8
3600.0 20 7

DataFrames implement a convenient plot that plots all columns against the DataFrame’s in-
dex. By default, it uses the matplotlib backend, so the plots can be modified further using
other Matplotlib functions. This is very useful to get a quite visual overview of the contained data.

# [Link](temperatures) # you can pass a dataframe to matplotlib and it␣


↪will plot all columns
[Link](); # but if you use pandas own wrapper it will also␣
↪add more useful labels

We can easily add new columns by assigning to a new key just like we would with a dictionary. If
we compute the new column from existing columns, they will inherit the same index, so the new
data integrates seamlessly into the existing DataFrame.

# we can access the series of the dataframe by name


print(temperatures["inside_temp"])
# we can easily compute new columns from existing ones and add them to␣
↪the dataframe
temperatures["difference"] = temperatures["inside_temp"] - temperatures[
↪"outside_temp"]
temperatures

time
0.0 17
600.0 18
1200.0 20
1800.0 21
2400.0 21
3000.0 19
3600.0 20
Name: inside_temp, dtype: int64

8.2. DataFrames as Collections of Series 105


Python for Engineers

inside_temp outside_temp difference


time
0.0 17 10 7
600.0 18 9 9
1200.0 20 11 9
1800.0 21 11 10
2400.0 21 10 11
3000.0 19 8 11
3600.0 20 7 13

There are (too) many ways to index into a DataFrame, here are some of the most common ones:

# directly indexing on the dataframe to get a column or a sub dataframe


temperatures[["outside_temp", "inside_temp"]] # the order matters

outside_temp inside_temp
time
0.0 10 17
600.0 9 18
1200.0 11 20
1800.0 11 21
2400.0 10 21
3000.0 8 19
3600.0 7 20

# use .loc to index by labels


[Link][1000:2000] # loc will use the index we defined (in this␣
↪case time in seconds)

inside_temp outside_temp difference


time
1200.0 20 11 9
1800.0 21 11 10

[Link][0:1800, "inside_temp"] # index row and column␣


↪simultaneously

time
0.0 17
600.0 18
1200.0 20
1800.0 21
Name: inside_temp, dtype: int64

# iloc will indexes by position like we know from most other array like␣
↪data structures
[Link][0:2, 1] # the first two rows, second column

time
0.0 10
600.0 9
Name: outside_temp, dtype: int64

8.2. DataFrames as Collections of Series 106


Python for Engineers

8.3 Convenience Functions

Pandas also implements a bunch of common statistical operations as member functions to Se-
riess and DataFrames.

# on Series
print(f"{temperatures["inside_temp"].mean() = }")
print(f"{temperatures["inside_temp"].sum() = }")
# ... and many more methods

temperatures["inside_temp"].mean() = np.float64(19.428571428571427)
temperatures["inside_temp"].sum() = np.int64(136)

# and on DataFrames
[Link]()

inside_temp 1.511858
outside_temp 1.511858
difference 1.914854
dtype: float64

# all (some) of them at once for the whole dataframe


[Link]()

inside_temp outside_temp difference


count 7.000000 7.000000 7.000000
mean 19.428571 9.428571 10.000000
std 1.511858 1.511858 1.914854
min 17.000000 7.000000 7.000000
25% 18.500000 8.500000 9.000000
50% 20.000000 10.000000 10.000000
75% 20.500000 10.500000 11.000000
max 21.000000 11.000000 13.000000

8.4 Indexing & Filtering

With Pandas, we get a lot of functionality that makes it easier to filter for the data you are interested
in. For example, let’s say we want to isolate all points in time where the inside temperature was
between 18 and 20 degrees:

# just like regular arithmetic operations, comparisons are also element-


↪wise
above_18 = 18 <= temperatures["inside_temp"]
# the comparison returns a boolean series indicating where the condition␣
↪is true
print(above_18)
below_20 = temperatures["inside_temp"] <= 20
# we can combine the boolean series with logical operators
indices = above_18 & below_20
# if we use a boolean series to index a dataframe we get only the rows␣
(continues on next page)

8.3. Convenience Functions 107


Python for Engineers

(continued from previous page)


↪where the condition is true
temperatures[indices]

time
0.0 False
600.0 True
1200.0 True
1800.0 True
2400.0 True
3000.0 True
3600.0 True
Name: inside_temp, dtype: bool

inside_temp outside_temp difference


time
600.0 18 9 9
1200.0 20 11 9
3000.0 19 8 11
3600.0 20 7 13

This is called boolean indexing and is also available for the underlying NumPy arrays.
One downside of filtering this way is that we have to repeat the name of the DataFrame for every
filter and indexing operation. This can make the code a bit ugly and harder to change. We can
get around this by using the query method.

filtered_boolean_indexing = temperatures[(18 <= temperatures["inside_temp


↪"]) & (temperatures["inside_temp"] <= 20)]
filtered_query = [Link]("18 <= inside_temp <= 20") # does␣
↪the same as the line above
filtered_boolean_indexing.equals(filtered_query)

True

8.5 Saving & Loading DataFrames

To do useful work on data, you will most likely need to load your data from somewhere and then
save processed results somewhere else. Pandas makes this very easy with its read_csv and
to_csv functions.

from [Link] import Pretty

temperatures.to_csv("[Link]")
# temperature_df.to_csv("[Link]", sep="\t") # the separator can␣
↪also be changed
Pretty(filename="[Link]")

time,inside_temp,outside_temp,difference
0.0,17,10,7
600.0,18,9,9
1200.0,20,11,9
(continues on next page)

8.5. Saving & Loading DataFrames 108


Python for Engineers

(continued from previous page)


1800.0,21,11,10
2400.0,21,10,11
3000.0,19,8,11
3600.0,20,7,13

When reloading the file, make sure to specify the correct index column, otherwise a new integer
index will be created, and the original index will just become another column.

reloaded_csv = pd.read_csv("[Link]", index_col="time")


reloaded_csv

inside_temp outside_temp difference


time
0.0 17 10 7
600.0 18 9 9
1200.0 20 11 9
1800.0 21 11 10
2400.0 21 10 11
3000.0 19 8 11
3600.0 20 7 13

While csv has its advantages like, ease of implementation, human readability and extreme preva-
lence and portability. It is not very space efficient or fast and has no standardized way of storing
metadata or preserving types. Luckily, Pandas also supports more advanced binary formats like
parquet, feather and hdf5 which are just as easy to use.

temperatures.to_parquet("[Link]")
reloaded_parquet = pd.read_parquet("[Link]") # no need to␣
↪specify the index column

8.6 Grouping & Aggregating

Finally, let’s take a look at how we can group and aggregate data in a DataFrame with a real
world example. The data we are working with is derived from this larger dataset.
Our contains the following columns for the months of April and May in 2015:
• datetime_str: The date and time of the observation as a string.
• temperature_C: The temperature in degrees Celsius.
• summary: A short description of the weather conditions.
We want to create four plots from this:
• The temperature over the entire time period.
• The daily minimum and maximum temperature and the difference between them over the
entire time period.
• A scatter plot of temperature over all hours of the day for May.
• A bar chart of the mean temperatures for the five most common weather conditions.

8.6. Grouping & Aggregating 109


Python for Engineers

weather = pd.read_csv("[Link]")

# process the datetime string into a datetime object and create helper␣
↪columns
weather["datetime"] = pd.to_datetime(weather["datetime_str"])
weather["date"] = weather["datetime"].[Link] # access individual␣
↪components of the datetime for all rows
weather["hour"] = weather["datetime"].[Link]
weather = weather.set_index("datetime").sort_index()

[Link](y=["temperature_degC"]);

To compute the daily minimum and maximum temperature, we will need to group the hourly data
by day and then compute the minimum and maximum temperature for each day. Since this is a
very common operation, Pandas has some very convenient methods for this. We start by grouping
the DataFrame by the auxiliary date column using the groupby method. It returns a GroupBy
object which is basically collection of tuples containing the groups name and a DataFrame con-
taining the data for that group.
for group, df in [Link]("date"):
print(f"{group = }")
display([Link](3)) # head(3) to only show the first 3 rows
break # only print the first group

group = [Link](2015, 4, 1)

datetime_str temperature_
↪degC \
datetime
2015-04-01 02:00:00+02:00 2015-04-01 02:00:00.000 +0200 9.
↪872222

2015-04-01 03:00:00+02:00 2015-04-01 03:00:00.000 +0200 8.


↪838889
2015-04-01 04:00:00+02:00 2015-04-01 04:00:00.000 +0200 8.
(continues on next page)

8.6. Grouping & Aggregating 110


Python for Engineers

(continued from previous page)


↪305556

summary date hour


datetime
2015-04-01 02:00:00+02:00 Overcast 2015-04-01 2
2015-04-01 03:00:00+02:00 Overcast 2015-04-01 3
2015-04-01 04:00:00+02:00 Overcast 2015-04-01 4

While we could work with this iterable ourselves now, we usually want to apply some common
aggregation function like mean, min, max, etc. to the groups. We do this by simply calling the
aggregations functions as member functions on the GroupBy object. It is also possible to perform
multiple aggregations at once by passing a list of aggregations to the agg method which is also a
member of the GroupBy object.

# only use the date and temperature columns because we only want
# to aggregate the temperature we need the date to group by
weather[["date", "temperature_degC"]].groupby("date").mean()

temperature_degC
date
2015-04-01 8.578535
2015-04-02 7.522454
2015-04-03 5.137037
2015-04-04 6.156481
2015-04-05 7.329630
... ...
2015-05-28 12.442130
2015-05-29 14.666898
2015-05-30 18.041667
2015-05-31 18.919213
2015-06-01 14.038889

[62 rows x 1 columns]

# it is also possible to select a column (or subset of columns) after␣


↪grouping
# we can compute multiple aggregates by passing a list of aggregation␣
↪function (names)

aggregates = [Link]("date")["temperature_degC"].agg(["min", "max


↪"])
aggregates["diff"] = aggregates["max"] - aggregates["min"]
display(aggregates)
[Link](rot=30)
[Link]("Temperature in °C");

min max diff


date
2015-04-01 4.250000 12.127778 7.877778
2015-04-02 1.738889 12.288889 10.550000
2015-04-03 1.627778 8.911111 7.283333
2015-04-04 -0.522222 12.127778 12.650000
2015-04-05 2.727778 11.066667 8.338889
... ... ... ...
2015-05-28 8.933333 16.138889 7.205556
(continues on next page)

8.6. Grouping & Aggregating 111


Python for Engineers

(continued from previous page)


2015-05-29 5.977778 22.222222 16.244444
2015-05-30 10.372222 24.972222 14.600000
2015-05-31 12.366667 24.727778 12.361111
2015-06-01 13.772222 14.305556 0.533333

[62 rows x 3 columns]

[Link]("[Link] == 5").plot(x="hour", y="temperature_degC


↪", kind="scatter");

# most operations in pandas can be chained together into a single␣


↪expression
most_common_summaries = weather["summary"].value_counts().head(5).keys()
(continues on next page)

8.6. Grouping & Aggregating 112


Python for Engineers

(continued from previous page)


[Link]("summary in @most_common_summaries").groupby("summary")[
↪"temperature_degC"].mean().plot(kind="bar", rot=0)
[Link]("Mean temperature in °C");

8.7 Summary

In this lecture we learned about:


• Pandas’ most important data structures the Series and the DataFrame.
• How to perform element wise operations on them.
• Filtering and aggregating data.
• Saving and loading DataFrames in different formats.
• More complex grouping and aggregation operations on datasets.
We really only had time for the absolute basics here, hopefully it was enough to give you a taste
and know where to start. Most of the things you would want to do with Pandas are probably already
implemented in some way and its certainly worth searching for these before hacking something
together yourself.
Some of the things we did not cover are:
• Working with unordered data.
• Transforming data with the apply method.
• Combining DataFrames with merge and concat.
• Reshaping data with pivot.
• And much more…

8.7. Summary 113


CHAPTER

NINE

ADVANCED PLOTTING

Open in JupyterLite: Solved Blank


One of the most important ways of communicating data and ideas is through visualization. And
while often times, a simple xy plot fulfills this purpose adequately, as the dimensionality of the
data increases, more advanced techniques are needed to keep the visualization concise and
informative.

import numpy as np
import [Link] as plt
from vscdark import set_style

set_style("vsclight")

<[Link] object>

9.1 Plots for Visualizing Distributions of Data

Up until now, we have only used line plots and scatter plots which are also really just line plots
without the lines. However, if we take a look at this Matplotlib cheatsheet we see that there are
many more types of plots that you might want to use. We are going to take a look at three more
of these, which are all associated with visualizing distributions in some way.
• [Link] to make histograms of univariate distributions of data.
• [Link] to compare multiple univariate distributions by characteristic values like the
median and quartiles.
• [Link] to compare multiple univariate distributions by approximate probability
density functions (kernel density estimates).
To compare different processes of producing 3 mm steel balls for ball bearings, we have measured
the diameter of 500 balls from each process.

[Link](42)
process_1 = [Link]([Link]((3, 3.01), (0.01, 0.02), (250,
↪ 2)))
process_2 = [Link](3.02, 0.005, 500)
process_3 = [Link]([[Link](3.0, 0.03, 490), [2.83, 2.85,
↪ 2.87, 2.89, 2.91, 2.93, 2.95, 2.97, 2.99]])

114
Python for Engineers

Histograms are a simple means of visualizing a distribution of values or approximating a proba-


bility density function from measured data. If the data is continuous or very sparse, we (or in this
case Matplotlib) bunch the values together in a set of bins that span the range of measured data.
We get an estimate of the distribution density function by dividing the number of measurements
in each bin by the total number of measurements and the width of the bin.
If you want something a bit smoother you might want to take a look at kernel density estimate plots
which are smoother but also potentially more misleading if you are not careful.

# there are a lot of options like specifying the number of bins,␣


↪normalizing to a probability density, etc.
[Link](process_1, bins=20, density=True)
[Link](False)
[Link]("Density")
[Link]("Diameter (mm)");

We can see that the distribution is slightly skewed to the left but otherwise mostly centered around
the desired diameter of 3 mm.
But how does it compare to the other two processes? To compare multiple distributions we can
plot them together in a single diagram using box plots or violin plots. Box plots and violin plots are
different ways of visualizing the same data: box plot show characteristic values like the median,
25% and 75% percentiles, and outliers, while violin plots also show a kernel density estimate of
the distribution.

fig, (axl, axr) = [Link](1, 2, figsize=(8, 5), sharey=True)


[Link]([process_1, process_2, process_3])
axl.set_xticklabels(["Process 1", "Process 2", "Process 3"])
axl.set_ylabel("Diameter (mm)")
[Link]([process_1, process_2, process_3])
axr.set_xticks([1, 2, 3], ["Process 1", "Process 2", "Process 3"]);

9.1. Plots for Visualizing Distributions of Data 115


Python for Engineers

Matplotlib only offers simple and sometimes fairly clunky to use versions of these plots. For more
advanced (and beautiful) statistical visualizations you should take a look at Seaborn. It is a sta-
tistical plotting library build on top of Matplotlib that abstracts away a lot of the more finicky details
of making these plots, offering a more high-level interface with good defaults.

9.2 Interactive Plots with Jupyter Widgets

Interactive plots can be very useful to give the reader an intuition of how changing a variable
can effect an otherwise static visualization. A user-changeable parameter can also serve as an
additional visualization dimension as an alternative or addition to 3D plotting. You can do a lot of
things with widgets, this article gives a good overview of most of the features.
The simplest way to use ipywidgets is to use interact. The basic idea is that we wrap the
thing we want to interact with in a function whose parameters we want to change interactively.
We then pass this function to interact, and it will create the necessary sliders, text boxes, and
other required input elements and automatically call the function with the current input elements’
values.

from ipywidgets import interact

def square(x: int = 2) -> int: # interact uses the function default to␣
↪initialize the slider
return x * x

interact(square, x=(-10, 10, 2)); # start, stop (, step)

interactive(children=(IntSlider(value=2, description='x', max=10, min=-


↪10, step=2), Output()), _dom_classes=('…

While the simple tuple syntax allows us to quickly get something up and running, we can get more
precise controls by passing the input widgets explicitly. And we are of course not limited to just
numbers or text as output, the function is allowed to produce or return any kind of output supported
in Jupyter notebooks including plots or HTML.

9.2. Interactive Plots with Jupyter Widgets 116


Python for Engineers

from ipywidgets import Dropdown, FloatSlider, Layout

def plot_trig(f: float, a: float, func):


x = [Link](0, 1, 100)
y = a * func(x * 2 * [Link] * f)
[Link](x, y)
[Link]([-2.1, 2.1])

interact(
plot_trig,
f=FloatSlider(value=1, min=0, max=4, step=0.1, layout=Layout(width=
↪"600px")),
a=FloatSlider(value=1, min=-2, max=2, step=0.1, layout=Layout(width=
↪"600px")),

func=Dropdown(options={"sin": [Link], "cos": [Link]}, value=[Link]),


);

interactive(children=(FloatSlider(value=1.0, description='f',␣
↪layout=Layout(width='600px'), max=4.0), FloatSli…

Using interact like this allows us to quickly create interactive visualizations. You can get a lot more
involved with different kinds of inputs and more complex user interfaces, but this is beyond the
scope of this lecture and for you to look up once you need it.

9.3 Matplotlib Interactive Backend

Another kind of interactivity you might want is the interactive backend of Matplotlib. Instead of
getting an image, you can also configure Matplotlib to return an interactive window which allows
you to zoom, pan and select data points. This is especially useful for large detailed datasets or
3D plots where you can rotate the plot to get a better understanding of the data.

# jupyter magic command to switch to the interactive backend


# uncomment for interactive mode
# %matplotlib widget
# [Link]() # disable automatic showing of plots

fig, ax = [Link]()

x = [Link](0, 2 * [Link], 100)


[Link](x, [Link](x));

# [Link]() # show figure explicitly

9.3. Matplotlib Interactive Backend 117


Python for Engineers

Matplotlib’s interactive backend also allows widgets to be added, but it is not as easy to use
as the ipywidgets, which should probably be your first choice for interactivity when working in
notebooks. You can find examples of this in the matplotlib documentation.
The interactive backend can be very useful, especially when exploring data, sadly it is not without
its caveats. Here are some downsides you might have to deal with:
• The interactive figures are kept open until they are explicitly closed. This can lead to perfor-
mance problems, especially when running cells multiple times.
• Since the interactive figures refer to live objects, they change when the object changes.
This can lead to unexpected behavior when using the same figure object in multiple cells
because later cells can change and even break the figures of earlier cells.
• Since the figures require the interpreter to be running, the outputs will be lost when reopening
the notebook and will need to be re-rendered. This way you can not simply view a notebook
without needing to re-execute the code.
• They are not very robust and may lead to errors on some platforms or with some configura-
tions.
So it might not be worth activating them when you don’t really need them.
If you want robust interactive plots that do not suffer from all of those problems, consider taking a
look at Plotly which outputs pure JavaScript and HTML that runs independently of the IPython
kernel.

9.3. Matplotlib Interactive Backend 118


Python for Engineers

9.4 3D Plotting

One application where the interactive backend really provides a lot of value is 3D plotting because
it allows rotating the plot to get a better view and 3D effect.

[Link]["[Link]"] = "azel" # disable roll

fig = [Link]()
ax = fig.add_subplot(projection="3d")

x = [Link](0, 4 * [Link])
s = 1j - 0.1
y = [Link]([Link](s * x))
z = [Link]([Link](s * x))

[Link](x, y, z)
ax.set_xlabel("x")
ax.set_ylabel("Real")
ax.set_zlabel("Imaginary")
ax.view_init(elev=20, azim=-50) # manually set the view for static output

# [Link]()

Instead of two outputs as a function of a single variable resulting in a curve in 3D space, we can
also plot a single function of two variables as a surface. To generate the grid of coordinates to
evaluate the function on, we can use [Link] which takes two (or more) 1D arrays and
returns two (or more) 2D arrays where each pair (or set) of values from the arrays corresponds to
a point on the grid.

9.4. 3D Plotting 119


Python for Engineers

[Link]([1, 2, 3, 4], [5, 6, 7])

(array([[1, 2, 3, 4],
[1, 2, 3, 4],
[1, 2, 3, 4]]),
array([[5, 5, 5, 5],
[6, 6, 6, 6],
[7, 7, 7, 7]]))

fig = [Link]()
ax = fig.add_subplot(projection="3d")

ls = [Link](-2, 2, 25)
x, y = [Link](ls, ls) # use meshgrid to evaluate a function on a␣
↪grid
z = [Link](-(x**2 + y**2))
print(f"{[Link] = }, {[Link] = }, {[Link] = }, {[Link] = }")

ax.plot_surface(x, y, z, cmap="viridis", alpha=0.5);


# color and transparency may make the plot easier to read
# [Link]()

[Link] = (25,), [Link] = (25, 25), [Link] = (25, 25), [Link] = (25,
↪ 25)

Keep in mind…
While 3D plots can look quite nice and have the potential to visualize things that 2D plots can’t,
they can quickly become confusing and are often hard to read. Always ask yourself whether you

9.4. 3D Plotting 120


Python for Engineers

really need a 3D plot or whether another kind of plot is better at getting your message across.
In most situations, the extra dimension can also be represented using color, size, or an extra
plot. Your number-one priority should always be readability. If your medium allows for it, use an
interactive plot or an animation so that the 3D effect is easier to grasp without perspective illusions.
Depending on your type of data, you can also draw lines to a reference plane or line.

# the same plot can also easily be created in 2D with color as the third␣
↪dimension
fig, ax = [Link]()
# cont = [Link](x, y, z, shading="gouraud") # full continuous␣
↪color
# cont = [Link](x, y, z, levels=20) # colored contour lines
cont = [Link](x, y, z, levels=20) # discrete levels make the␣
↪contours easier to see
[Link](cont)
[Link](False)
# [Link]()

9.5 Animations

Another dimension we can make use of is time itself by animating our plots. For this, matplotlib
provides the FuncAnimation class.

from [Link] import FuncAnimation

fig, ax = [Link]()

x = [Link](0, 2 * [Link])
line = [Link](x, [Link](x))[0] # save the line handle from the plot call

N = 100
(continues on next page)

9.5. Animations 121


Python for Engineers

(continued from previous page)

def animate(frame_num):
line.set_ydata([Link](x - frame_num / N * 2 * [Link]))
return (line,) # return all the objects that we changed

anim = FuncAnimation(fig, animate, N, interval=20)


# [Link]()

Displaying these animations live obviously requires the interactive backend. But like other uses
of the interactive backend, this can be a bit unreliable. One work around for this is to save the
file to a video or GIF and display that instead. This also has the advantages that the animations
persists when reopening the notebook and can be shared more easily. However, it can be a bit
slow to render.

from [Link] import HTML, Video

# either directly create html with embedded video


HTML(anim.to_html5_video())

# or save the animation to a file and display it


# [Link]("sine.mp4")
# Video(filename="./sine.mp4", embed=True)

<[Link] object>

FuncAnimation is also not without its downsides and technicalities. Depending on your work-
flow, it might also make sense to generate and save the individual frames one by one and then
use something like FFmpeg to combine them into a video.

9.5. Animations 122


Python for Engineers

9.6 Summary

In this lecture we learned about:


• Creating interactive plots and other interactive outputs using ipywidgets’ interact
function.
• Using Matplotlib’s interactive backend to explore your data.
• Creating 3D plots using [Link] and Matplotlib’s 3D plotting capabilities.
• How to create animations and embed them in your notebook.

9.6. Summary 123


CHAPTER

TEN

DIFFERENTIAL EQUATIONS & SYSTEMS

Open in JupyterLite: Solved Blank


Differential equations are a very natural way of describing our physical reality. But while they are
often fairly easy to formulate, they are often very hard or even impossible to solve analytically.
Luckily, we live in a world of plentiful computational resources and Python makes it very easy to
exploit these resources to numerically solve differential equations.
Differential equations can be categorized into two main types: ordinary differential equations
(ODEs) and partial differential equations (PDEs). ODEs only involve functions and derivatives
with respect to one variable like a harmonic oscillator:

𝑑2 𝑥(𝑡)
+ 𝜔2 𝑥(𝑡) = 0
𝑑𝑡2
Here 𝑥(𝑡) is only a function of time and the derivative is also only with respect to time 𝑡. In these
2
cases, we typically use the Newton notation for the time derivatives, e.g., 𝑥̈ = 𝑑𝑑𝑡2𝑥 for second
derivative.
Meanwhile, PDEs involve functions and derivatives with respect to multiple variables, such as the
heat equation:

𝜕𝑢(𝑡, 𝑥) 𝜕 2 𝑢(𝑡, 𝑥)
=𝛼
𝜕𝑡 𝜕𝑥2
Here 𝑢(𝑡, 𝑥) is a function of both time 𝑡 and space 𝑥 and the derivatives are partial derivatives.
In this chapter, we will only cover ODEs, but the methods we will learn can be extended to PDEs
as well. This can be done by approximating the continuous space with a discrete grid which yields
a system of ODEs.

import numpy as np
import [Link] as plt

[Link]("[Link]")

124
Python for Engineers

10.1 First Order Differential Equations

Let’s consider this simple non-linear differential equation:

𝑦(𝑡)
̇ = −𝑦(𝑡) ⋅ 𝑡

with initial the condition 𝑦(0) = 1.


This is an initial value problem (IVP), meaning we know the value for all degrees of freedom at a
single point in the functions domain, in this case time. This is opposed to a boundary value problem
(BVP) where the constraints are given at different points in the functions’ domain. Depending
on the so-called “boundary conditions”, a boundary value problem can have multiple or even no
solutions, while an initial value problems always have one unique solution. We will only consider
initial value problems in this chapter.
Since it is an IVP, we will use the appropriately named solve_ivp function from the scipy.
integrate module to solve it.
There are three mandatory arguments for solve_ivp(fun, t_span, y0, ...):
• fun is how we are supposed to provide the differential equation to the solver.
• t_span is the time span that we want to solve the differential equation over.
• y0 is the initial condition.
The first argument is probably the most interesting one, it asked us to supply a function that satis-
fies 𝑑𝑦
𝑑𝑡 = 𝑓(𝑡, 𝑦). What is this function for us? It is simply the right side of the differential equation:

!
𝑦(𝑡)
̇ = −𝑦(𝑡) ⋅ 𝑡 = 𝑓(𝑡, 𝑦)

⇒ 𝑓(𝑡, 𝑦) = −𝑦 ⋅ 𝑡
Another very useful optional argument is t_eval, which can be used to specify at which points
in time the solution should be determined. Without this, the solver will choose the points in time
itself, attempting to keep the error below a certain threshold.

from [Link] import solve_ivp

delta_t = 0.1
ts = [Link](0, 3, delta_t)
interval = (min(ts), max(ts))
y0 = [1]

# let the solver decide the evaluation points


sol = solve_ivp(lambda t, y: -y * t, interval, y0)
# limit the step size
sol_max_step = solve_ivp(lambda t, y: -y * t, interval, y0, max_
↪step=delta_t)

# explicitly specify the (returned) evaluation points


sol_t_eval = solve_ivp(lambda t, y: -y * t, interval, y0, t_eval=ts)
[Link](sol.t, sol.y[0], label="solve_ivp (default)")
[Link](ts, [Link](-(ts**2) / 2), label="analytical solution")
[Link](sol_max_step.t, sol_max_step.y[0], label="solve_ivp (max_step)",␣
↪linestyle="--")

[Link](sol_t_eval.t, sol_t_eval.y[0], label="solve_ivp (t_eval)",␣


↪linestyle=":")
[Link]();

10.1. First Order Differential Equations 125


Python for Engineers

There are many more options to solve_ivp. These are including but not limited to the solver to
use, tolerance settings, what exactly to output and solver specific additional information about the
problem to speed up the process. If your results are incorrect or too slow, you should probably
take a look at these settings. However, the default settings should be sufficient for most cases.

10.2 Higher Order Differential Equations

Now let’s look at something slightly more complicated, a harmonic oscillator. A simple driven
harmonic oscillator can be described by the following equation:

𝑚𝑦(𝑡)
̈ + 𝑑𝑦(𝑡)
̇ + 𝑘𝑦(𝑡) = 𝑢(𝑡)

𝑑 𝑘 1
⇒ 𝑦(𝑡)
̈ =− 𝑦(𝑡)
̇ − 𝑦(𝑡) + 𝑢(𝑡)
𝑚 𝑚 𝑚
where 𝑚, 𝑑, and 𝑘 are the mass, damping coefficient, and spring constant respectively. 𝑦(𝑡) is the
displacement of the mass at time 𝑡, and 𝑢(𝑡) is an external force applied to the mass.
If we want to plug this into solve_ivp, we quickly notice a problem. The function requires us
to provide a function that returns the first derivative of the variable to integrate, but there is a
second derivative on the left side. While solve_ivp can actually only solve equations with first
order time derivatives, it can also solve entire systems of first order differential equations. This is
not a big problem however, since any 𝑛-th order ODE can be reformulated into a system of 𝑛 first
order ODEs by introducing new variables for the lower order derivatives. In this case we introduce
the velocity 𝑣(𝑡) = 𝑦(𝑡)
̇ as a new variable so that 𝑦(𝑡) ̈ = 𝑣(𝑡)
̇ is only a first order derivative of the
state.
𝑦(𝑡)
̇ 𝑦(𝑡)
̇ 𝑦(𝑡)
̇ =[
x(𝑡) ]=[1 𝑑 𝑘 ] = f (𝑡, x = [ ])
𝑣(𝑡)
̇ 𝑚 𝑢(𝑡) − 𝑚 ̇
𝑦(𝑡) − 𝑚 𝑦(𝑡)
𝑣(𝑡)

For now let’s consider the autonomous case 𝑢(𝑡) = 0 and see how the system evolves from the
initial conditions 𝑦(0) = 1 and 𝑦(0)
̇ = 0.

10.2. Higher Order Differential Equations 126


Python for Engineers

k = m = 1
d = 0.3

def dydt_oscillator(t: float, ys: [Link]) -> [Link]:


y, y_dot = ys # y is now a vector
return [y_dot, -k / m * y - d / m * y_dot]

interval = (0, 10)


dt = 0.1

sol = solve_ivp(dydt_oscillator, interval, [1, 0], max_step=dt)


[Link](sol.t, sol.y[0], label="y")
[Link](sol.t, sol.y[1], label="v")
[Link]();

10.3 Differential Equations with Inputs

Now let’s consider the case where we have a harmonic excitation 𝑢(𝑡) = sin(𝜔𝑡). We can imple-
ment this excitation as a time dependent term in the dydt function.
from typing import Callable

dt = 0.1
interval = (0, 50)

# naive way, baked into the function


# f = 1
# def dydt(t: float, ys: [Link]) -> [Link]:
# y, y_dot = ys # y is now a vector
# return [y_dot, [Link](2 * [Link] * f * t) - k / m * y - d / m * y_
↪dot]

# sol = solve_ivp(dydt, interval, [0, 0], max_step=dt)


# [Link](sol.t, sol.y[0], label="y")
(continues on next page)

10.3. Differential Equations with Inputs 127


Python for Engineers

(continued from previous page)

# more flexible, pass the function as an argument to dydt


def dydt_oscillator_input(t: float, ys: [Link], u: Callable[[float],␣
↪float]) -> [Link]:
y, y_dot = ys # y is now a vector
return [y_dot, u(t) - k / m * y - d / m * y_dot]

for w in [0.5, 1, 2]:


u = lambda t, w=w: [Link](w * t)
# additional arguments to `dydt` (besides t and y) can be passed␣
↪using `args`
sol = solve_ivp(dydt_oscillator_input, interval, [0, 0], args=(u,),␣
↪max_step=dt)
[Link](sol.t, sol.y[0], label=f"y (u=sin({w}t))")

[Link]();

10.4 Linear Time-Invariant Systems

A lot of the systems we encounter are (approximately) linear and time-invariant systems. This
means we can make use of all the tools that you learned in your control theory and signals and
systems class. To make use of these tools in Python, we will use the Python Control Systems
Library often called control for short. Many of these features are also implement in scipy.
signal which also contains a lot of useful functions for signal processing, but we will focus on
control here.

10.4. Linear Time-Invariant Systems 128


Python for Engineers

10.4.1 State Space Representation

The driven harmonic oscillator from the previous section is a linear time-invariant system because
the system function f(𝑡, x) is linear in x(𝑡) and 𝑢(𝑡) and otherwise independent of 𝑡. When this is
the case, it can be written as a matrix vector expression:

̇ = f(𝑡, x(𝑡)) = Ax(𝑡) + B𝑢(𝑡)


x(𝑡)

we also introduce the concept of an output 𝑦(𝑡) which may not necessarily include all states of the
system:

𝑦(𝑡) = Cx(𝑡) + 𝐷𝑢(𝑡)

For the driven harmonic oscillator we get:

0 1 0
̇ =[
x(𝑡) 𝑘 𝑑 ] x(𝑡) + [ 1 ] 𝑢(𝑡)
−𝑚 −𝑚 𝑚
𝑦(𝑡) = [1 0] x(𝑡) + 0 ⋅ 𝑢(𝑡)

To work with this in Python Control, we use the StateSpace class:

import control as ct

A = [Link]([[0, 1], [-k / m, -d / m]])


B = [Link]([[0], [1 / m]])
C = [Link]([[1, 0]])
D = [Link]([[0]])

harmonic_oscillator_sys = [Link](A, B, C, D)
harmonic_oscillator_sys

<StateSpace sys[0]: [‘u[0]’] -> [‘y[0]’]>

0 1 0
⎡ −1 −0.3 1 ⎤
⎢ ⎥
⎣ 1 0 0 ⎦

The poles of the system, which are equal to the eigenvalues of the matrix A, show us that the
system is stable with our choice of parameters.

print(f"{harmonic_oscillator_sys.poles() = }")
print(f"{[Link](A) = }")
print("stable?", [Link]([Link](harmonic_oscillator_sys.poles()) < 0))

harmonic_oscillator_sys.poles() = array([-0.15+0.988686j, -0.15-0.


↪988686j])
[Link](A) = array([-0.15+0.988686j, -0.15-0.988686j])
stable? True

We can excite the system with an arbitrary sequence of inputs. Let’s use the same harmonic
excitation to make sure we get the same results.

10.4. Linear Time-Invariant Systems 129


Python for Engineers

for w in [0.5, 1, 2]:


ts = [Link](*interval, dt)
us = [Link](w * ts)
res = ct.forced_response(harmonic_oscillator_sys, T=ts, U=us)
[Link]([Link], [Link], label=rf"y (u=sin({w}t))")
[Link]();

To get the resulting amplitude and phase lag for a whole frequency range, we can use the fre-
quency_response function. The result is basically a numerical representation of the system’s
transfer function.

# the frequency range is computed automatically based on the systems␣


↪poles and zeros
resp = harmonic_oscillator_sys.frequency_response()
print(f"{type(resp) = }")
[Link]();
# sys.bode_plot() # or use the bode_plot method which also calculates the␣
↪frequency response

type(resp) = <class '[Link]'>

10.4. Linear Time-Invariant Systems 130


Python for Engineers

Two other signatures of a system are its impulse and step response. They can be computed using
the impulse_response and step_response functions respectively. They return TimeRe-
sponseData objects which contains the time points and states/outputs of the system at those
time points. They also implement convenient plot methods to visualize the results.

step_response = harmonic_oscillator_sys.step_response()
print(f"{type(step_response) = }")
[Link](step_response.time, step_response.outputs, label="step response")
# sys.step_response().plot(label="step response")
harmonic_oscillator_sys.impulse_response().plot(label="impulse response",␣
↪color="C1");

type(step_response) = <class '[Link]'>

10.4. Linear Time-Invariant Systems 131


Python for Engineers

10.4.2 Transfer Function Representation

Another common representation for LTI systems are transfer functions. The control systems
library implements these using the TransferFunction class. We can create a transfer function
either by specifying the numerator and denominator polynomials or by specifying the zeros, poles
and a gain k using zpk.

[Link]([2], [1, 2, 1])

<TransferFunction sys[2]: [‘u[0]’] -> [‘y[0]’]>


2
𝑠2 + 2𝑠 + 1
[Link]([], [-1, -1], 2)

<TransferFunction sys[3]: [‘u[0]’] -> [‘y[0]’]>


2
𝑠2 + 2𝑠 + 1

There are also member functions to convert a StateSpace model to a TransferFunction


and vice versa.

harmonic_oscillator_sys.to_tf()

<TransferFunction sys[4]: [‘u[0]’] -> [‘y[0]’]>


1
𝑠2 + 0.3𝑠 + 1

10.4.3 Other Operations on LTI Systems

Control systems are usually made up of multiple LTI systems. Systems can be combined by
putting them in series which corresponds to multiplication of the transfer function or in parallel
which corresponds to addition of the transfer function.

plant = [Link]([2], [1, 1])


sensor = [Link]([1], [1 / 5, 1])

print("series")
display(plant * sensor)

print("parallel")
display(plant + sensor)

series

<TransferFunction sys[7]: [‘u[0]’] -> [‘y[0]’]>


2
0.2𝑠2 + 1.2𝑠 + 1

10.4. Linear Time-Invariant Systems 132


Python for Engineers

parallel

<TransferFunction sys[8]: [‘u[0]’] -> [‘y[0]’]>


1.4𝑠 + 3
0.2𝑠2 + 1.2𝑠 + 1

We can also add a feedback path using the feedback function or realize any other arbitrary
topology using the interconnect function.

integral_controller = [Link]([], [0], 1) # one pole at 0 is an integrator


open_loop = plant
closed_loop = (integral_controller * plant).feedback(sensor, sign=-1)
ct.step_response([open_loop, closed_loop]).plot(label=["open loop",
↪"closed loop"]);

The integral controller makes sure that the output matches the input eventually. In other words, it
ensures a DC-gain of 1 and thereby a steady state error of 0.
Some of these features are not limited to linear systems, they also work on potentially nonlinear
systems which are represented by the generic InputOutputSystem class. Non-linear systems,
discrete time systems and the much more is beyond the scope of this showcase, however.

10.4. Linear Time-Invariant Systems 133


Python for Engineers

There is a lot more you can do with the control library, but this is for you to look up as you need
it.

10.5 Phase Portraits

While step responses and frequency responses are a very useful tool for understanding a system’s
behavior, they are limited to LTI systems. Another insightful way to visualize system behavior that
also works for non-linear systems is the phase portrait.
As we saw earlier when solving non-linear differential equations with solve_ivp, the au-
tonomous behavior of a system can be described by the function x(𝑡) ̇ = f(𝑡, x(𝑡)). If the system is
time-invariant, it basically tells us the direction (and magnitude) that the system will change with
̇ for a given state x(𝑡). If the system has two states, this mapping can be visualized as a vector
x(𝑡)
field in the plane.
There are two kinds of matplotlib plots that are suited for visualizing this kind of data:
• quiver allows us to show both the magnitude and direction of the flow, but it is hard to
follow the actual trajectories in phase space.
• streamplot visualizes the trajectories or “streamlines” in phase space, we can still use
color or line width to show the magnitude of the flow.

ls = [Link](-1, 1, 20)
xs, vs = [Link](ls, ls)
# if you write your function right, they work with any dimension of numpy␣
↪arrays
# so we can easily evaluate the function on the entire grid
x_dots, v_dots = dydt_oscillator(0, [xs, vs])
magnitude = [Link]([x_dots, v_dots], axis=0)

fig, (axl, axr) = [Link](1, 2, figsize=(10, 5))


[Link](xs, vs, x_dots, v_dots, magnitude, pivot="mid")
[Link](False)
axl.set_xlabel("position")
axl.set_ylabel("velocity")
[Link](xs, vs, x_dots, v_dots, color=magnitude)
axr.set_xlabel("position")
[Link](False)
[Link]("Phase Portrait of Harmonic Oscillator");

10.5. Phase Portraits 134


Python for Engineers

We can see how a positive position will lead to a negative change in velocity, and a positive velocity
will lead to a positive change in position so that the trajectories in phase space are spirals around
the origin. The streamlines from streamplot also enable us to see that the trajectories spiral
inwards towards the origin, which is the stable equilibrium point of the system.
While it is good to be aware of these two kinds of plots, you should know that this functionality
and much more are already implemented in control.phase_plane_plot. In addition to the
streamlines or vectorfield, it can also compute and plot some additional features that are especially
interesting for non-linear systems.
A classic example of a simple non-linear differential equation is a pendulum:
̈ = − 𝑔 sin(𝜃(𝑡))
𝜃(𝑡)
𝑙
or in state space form:
̇
𝜃(𝑡) 𝜔(𝑡)
̇ =[
x(𝑡) ]=[ 𝑔 ]
𝜔(𝑡)
̇ − 𝑙 sin(𝜃(𝑡))
where 𝜃(𝑡) is the angle of the pendulum and 𝜔(𝑡) is the angular velocity.
g = l = 1

def dydt_pendulum(t: float, ys: [Link]) -> [Link]:


theta, omega = ys
return [omega, -g / l * [Link](theta)]

[Link](figsize=(8, 6))
ct.phase_plane_plot(dydt_pendulum, pointdata=[-4, 4, -3, 3], plot_
↪streamplot={"vary_linewidth": True})
[Link](r"angle $\theta$")
[Link](r"angular velocity $\omega$")
# ct.phase_plane_plot(harmonic_oscillator_sys) # also works with␣
↪StateSpace and other InputOutputSystem objects
[Link]().set_aspect("equal")
[Link](False)

10.5. Phase Portraits 135


Python for Engineers

Not only do we get the streamlines, phase_plane_plot also finds the equilibrium points and
plots separatrices around equilibrium points with real valued eigenvalues. As we can see, the
pendulum has stable (if it is dampened) equilibrium points at (𝜃 = 0 mod 2𝜋, 𝜔 = 0) and saddle
points at (𝜃 = 𝜋 mod 2𝜋, 𝜔 = 0) which separate trajectories that swing back and forth from those
that have enough energy swing around.

10.6 Summary

In this lecture we learned about:


• Solving differential equations using SciPy’s solve_ivp.
• Reformulating higher order differential equations as systems of first order differential equa-
tions to solve them numerically.
• Solving time dependent or inhomogeneous differential equations.
• Using the Python Control Systems Library to work with and analyze linear time-invariant
systems.

10.6. Summary 136


CHAPTER

ELEVEN

ROOT FINDING & OPTIMIZATION

Open in JupyterLite: Solved Blank


Many engineering related problems require us to solve or optimize complex equations or expres-
sion. Often times, these problems have no closed form solution or are otherwise too complex to be
solved analytically (even assisted by computer algebra systems like SymPy). In these cases we
usually resort to numerical methods to find approximate solutions to specific instances (i.e. with
given values) of these problems. In this lecture we will take a look at SciPy’s optimize module,
which provides a wide range of functions for root finding, optimization, and curve fitting.

import numpy as np
import [Link] as plt

[Link]("[Link]")

11.1 Solving Equations (Root Finding)

SciPy and most other numerical packages in other languages usually don’t directly provide func-
tion to solve equations, they only provide functions to find the roots (or zeros) of a function, i.e.
finding the values of 𝑥 for which 𝑓(𝑥) = 0. We can reformulate any equation solving problem
lhs(𝑥) = rhs(𝑥) into a root finding problem by subtracting the right-hand side from the left-hand
side:

𝑓(𝑥) = lhs(𝑥) − rhs(𝑥) = 0

Finding the roots or zeros of this function is equivalent to finding the solution to the original equa-
tion.

11.1.1 Scalar Problems

One of the simplest equations with no closed form solution that appears in real world engineering
problems is the following:

𝑒𝑥 = −𝑥

Equations of this form appear when solving for the current in a simple circuit with a diode attached
to a real voltage source.
𝑈
𝐼0 ⋅ (𝑒 𝑈𝑇 − 1) = (𝑈0 − 𝑈 )/𝑅𝑖

137
Python for Engineers

There is no way to express the solution to this equation in terms of elementary functions. However,
we can still find an arbitrarily precise* numerical solution for a given set of parameters by making
use of numerical root finding methods.
Reformulated to a root finding problem we get:
!
𝑓(𝑥) = 𝑒𝑥 + 𝑥 = 0

Let’s see how we can solve this in Python using functions from [Link]. The submod-
ule contains a few functions that we could use, but the go-to candidate is usually root_scalar,
which provides a unified interface to scalar root finding algorithms.

from scipy import optimize

f = lambda x: [Link](x) + x
res_init = optimize.root_scalar(f, x0=0) # we need to provide an initial␣
↪guess

res_bracket = optimize.root_scalar(f, bracket=[-1, 1]) # or a range


print("res_init:\n", res_init, sep="")
print("res_bracket:\n", res_bracket, sep="")

xs = [Link](-2, 2, 100)
[Link](xs, [Link](xs), label="lhs(x) = exp(x)")
[Link](xs, -xs, label="rhs(x) = -x")
[Link](res_init.root, [Link](res_init.root), color="red", label=
↪"Intersection", zorder=10)
[Link]();

res_init:
converged: True
flag: converged
function_calls: 10
iterations: 5
root: -0.567143290409784
method: newton
res_bracket:
converged: True
flag: converged
function_calls: 7
iterations: 6
root: -0.5671432904097838
method: brentq

11.1. Solving Equations (Root Finding) 138


Python for Engineers

11.1.2 Limitations

A critical caveat of these numerical methods is that most of them generally only find a root, not
all roots. The latter is a much harder problem and usually requires more knowledge about the
function you are trying to solve.
There are tricks that we can use to find multiple roots, but they usually require additional informa-
tion and have their own caveats and limitations.

xs = [Link](-2, 2, 100)
f_2 = lambda x: x**2 - 1

[Link](xs, f_2(xs), label="$f_2(x) = x^2 - 1$")

# only one solution per call


# one_sol = optimize.root_scalar(f_2, x0=0)
# [Link](one_sol.root, f_2(one_sol.root), color="red", label="One␣
↪solution", zorder=10)

# find many solutions by dividing the interval of interest into many␣


↪search regions ...
solutions = []
search_grid = [Link](-2, 2, 10)
for left, right in zip(search_grid[:-1], search_grid[1:], strict=True):
if f_2(left) * f_2(right) > 0:
continue
sol = optimize.root_scalar(f_2, bracket=[left, right])
[Link]([Link])

[Link](solutions, np.zeros_like(solutions), color="red", label=


↪"Solutions", zorder=10)
[Link]();

11.1. Solving Equations (Root Finding) 139


Python for Engineers

Dividing the interval of interest into smaller search regions does give us multiple roots, but there
are still lots of caveats. We can only find one root per segment, so the search grid needs to be fine
enough, another critical parameter that we need to know beforehand. Using the sign to determine
whether a root exists in an interval can fail the function changes sign an even number of times in
the interval. Lastly, the process can quickly become very computationally demanding, especially
in higher dimensions.

11.1.3 Multiple Unknowns and Equations

We can also solve systems of equations in a similar way. Just like we can reformulate a single
equation as a root finding problem 𝑓(𝑥) = 0, we can reformulate a system of equations as a
vectorial root finding problem f(x) = 0.
For example, let’s consider the following system of equations:

𝑥1 + 𝑥2 = 𝑥21 + 1
𝑥1 ⋅ 𝑥2 = 6

We can reformulate this in terms of the vector x = [𝑥1 , 𝑥2 ]𝑇 and the vectorial function f(x) =
[𝑓1 (x), 𝑓2 (x)]𝑇 :

𝑓1 (x) 𝑥 + 𝑥2 − 𝑥21 − 1
f(x) = [ ]=[ 1 ]
𝑓2 (x) 𝑥1 ⋅ 𝑥2 − 6

def f_vec(x: tuple[float, float]) -> tuple[float, float]:


x1, x2 = x
return (
x1 + x2 - x1**2 - 1,
x1 * x2 - 6,
)

res_vec = [Link](f_vec, x0=(0, 0))


print(res_vec)
(continues on next page)

11.1. Solving Equations (Root Finding) 140


Python for Engineers

(continued from previous page)


x1, x2 = res_vec.x
print(f"{x1 + x2} = {x1**2 + 1}")
print(f"{x1 * x2} = {6}")

message: The solution converged.


success: True
status: 1
fun: [-9.770e-14 3.242e-13]
x: [ 2.000e+00 3.000e+00]
method: hybr
nfev: 32
fjac: [[ 7.562e-01 -6.544e-01]
[-6.544e-01 -7.562e-01]]
r: [-3.750e+00 -9.228e-01 -2.402e+00]
qtf: [-7.685e-10 -4.860e-10]
5.000000000000133 = 5.000000000000231
6.000000000000324 = 6

11.2 Optimization

Optimization is very much related to root finding, since, as we know from calculus, a minimum or
maximum is a root of the derivative of the respective function (assuming the function is differen-
tiable). However unlike zeroes, not all stationary points are equal, some might be global minima
others might be local maxima and yet others might no minimum or maximum at all.
Let’s start with a simple example that we can also solve by just looking at the function.
We want to find the minimum of 𝑓3 (𝑥) = (𝑥 − 1)2
For this, we will use minimize and just like with the local root finding methods, we also need to
provide an initial guess.

def f_3(x: float) -> float:


return (x - 1) ** 2

[Link](f_3, x0=0)

message: Optimization terminated successfully.


success: True
status: 0
fun: 5.5507662238258444e-17
x: [ 1.000e+00]
nit: 2
jac: [ 4.682e-13]
hess_inv: [[ 5.000e-01]]
nfev: 6
njev: 3

It finds the obvious minimum at 𝑥 = 1.


Now let’s try something a little more tricky, we want to find the maximum of 𝑓4 (𝑥) = cos(2𝑥) −
1 2
4𝑥 + 𝑥

11.2. Optimization 141


Python for Engineers

The immediate problem that we encounter it that there is no maximize function. Instead, we
have to find the maximum by searching for a minimum of the negated function 𝑓4,neg = −𝑓4 (𝑥).

def f_4(x: float) -> float:


return [Link](2 * x) - x**2 / 4 + x

res_loc = [Link](lambda x: -f_4(x), x0=0)


res_loc

message: Optimization terminated successfully.


success: True
status: 0
fun: -1.1128285242030254
x: [ 2.293e-01]
nit: 4
jac: [-1.490e-08]
hess_inv: [[ 2.446e-01]]
nfev: 12
njev: 6

The optimization terminated successfully, but did we actually find the maximum? Let’s plot the
function in our region of interest and see what it looks like.

xs = [Link](-2, 5)
[Link](xs, f_4(xs))
[Link](res_loc.x, f_4(res_loc.x), color="red", zorder=10);

It seems like we only found a local maximum, which is a fundamental downside to local optimiza-
tion. Just like numerical root finding only finds a single root, numerical optimization only finds a
single optimum. But unlike roots, an optimum might only be local. If we had chosen a different ini-
tial guess we might have found the true optimum but choosing the right initial guess either requires
additional knowledge about the function.
Global optimization attempts to solve this problem by systematically searching the region of inter-
est instead of just iteratively improving an initial guess.

11.2. Optimization 142


Python for Engineers

One global optimizer provided by [Link] is differential_evolution. Like im-


plied in the previous sentence, we don’t supply an initial guess but rather a region that we want to
search.

res_glob = optimize.differential_evolution(lambda x: -f_4(x), bounds=[(-2,


↪ 5)])
res_glob

message:
Optimization terminated successfully.
success:
True
fun:
-1.7105680653395763
x:
[ 3.014e+00]
nit:
4
nfev:
81
population:
[[ 3.015e+00]
[ 3.022e+00]
...
[ 2.934e+00]
[ 3.020e+00]]
population_energies: [-1.711e+00 -1.710e+00 ... -1.697e+00 -1.710e+00]
jac: [ 2.220e-08]

It seems to have found the global optimum at 𝑥 ≈ 3, but it is important to know that this is not
guaranteed. Unlike the name might lead you to believe, global optimizers don’t necessarily find
the global optimum, they just search “globally”.
As always, there is much more to know here, just like the root finding functions, the optimization
functions also have a lot of options like:
• method: Choose a different optimization algorithm.
• bounds: Limit the region of valid solutions.
• jac, hess: Provide a way to compute the Jacobian or Hessian of the function to speed up
optimization. These can be automatically computed using libraries like jax.
• constraints: Account for constraints in the optimization problem.
• Various optimizer specific settings that impact performance and accuracy…

11.3 Curve Fitting

A very common kind of optimization problem that you likely already have encountered is curve
fitting. We usually have some input output pairs (𝑥𝑛 , 𝑦𝑛 ), and we want to find some parameters
𝑝1 , 𝑝2 , … such that a function 𝑓(𝑥𝑛 , 𝑝1 , 𝑝2 , …) ≈ 𝑦𝑛 fits our observations as well as possible. But
since “as best as possible” is not really something we can give to minimize we need to define
a loss function (also called cost or objective function) that quantitatively describes how good (or
technically speaking bad) our fit is.
A very common loss function for these kinds of problems is the mean-squared error (MSE) where
we just sum up the squared errors for all observations:
𝑁
𝐿(𝑝1 , 𝑝2 , …) = ∑(𝑓(𝑥𝑛 , 𝑝1 , 𝑝2 , …) − 𝑦𝑛 )2
𝑛=1

11.3. Curve Fitting 143


Python for Engineers

Technically this is not the mean squared error but the sum of squared errors, but this does not
change the result of the optimization.

# mock measurement data


[Link](0)
a0 = 2
gamma = 0.2
freq = 0.5
sigma = 0.1
ts = [Link](0, 10, 100)
ys = a0 * [Link](2 * [Link] * freq * ts + 0.5) * [Link](-gamma * ts) + np.
↪[Link](0, sigma, len(ts))

Say we measured the impulse response of a harmonic oscillator. We know it should be a decaying
sinusoidal function like this:

𝑦(𝑡) = 𝑎0 ⋅ sin(2𝜋𝑓𝑡 + 𝜑0 ) ⋅ 𝑒−𝛾𝑡

Let’s fit the function to our data to determine the dampened resonance frequency 𝑓 and the damp-
ing factor 𝛾. We will also need to fit the amplitude 𝑎0 .

[Link](ts, ys, "o-", label="Measurements")

def fit_func(t: float, a0: float, freq: float, gamma: float, phi: float) -
↪> float:
return a0 * [Link](2 * [Link] * freq * t + phi) * [Link](-gamma * t)

def loss_func(params: [Link]) -> float:


return [Link]((fit_func(ts, *params) - ys) ** 2)

params = [Link](loss_func, x0=[1, 1, 1, 0]).x


a0, freq, gamma, phi = params
[Link](ts, fit_func(ts, a0, freq, gamma, phi), label="Fit")
[Link]()
print(f"{a0 = }, {freq = }, {gamma = }, {phi = }")

a0 = np.float64(2.0657636218450524), freq = np.float64(0.


↪4994305663301355), gamma = np.float64(0.21281031849873955), phi = np.
↪float64(0.5026055376198153)

11.3. Curve Fitting 144


Python for Engineers

While this was a good exercise in formulating an objective functions, it is an unnecessarily com-
plicated approach. Since this is such a common task, [Link] already provides a
function curve_fit that formulates and solved the optimization problem for you.

params_curve_fit, cov = optimize.curve_fit(fit_func, ts, ys, p0=[1, 1, 1,␣


↪0])

# the results are pretty much identical


print(params)
print(params_curve_fit)

[2.06576362 0.49943057 0.21281032 0.50260554]


[2.06576378 0.49943058 0.21281035 0.5026053 ]

11.4 Optimal Control

On domain where we can make use of these numerical methods is optimal control. As the name
suggests, optimal control concerns itself with optimally controlling systems. This can be either
offline by finding optimal rules to control a system like optimal gains for a PID controller, or online by
finding optimal control inputs to make a system follow a desired trajectory as closely as possible.

11.4.1 Offline Optimization

For example, we may want to optimize the gains of a feedback controller so that it quickly reaches
its setpoint with minimal ringing. Let’s use the harmonic oscillator that we introduced last lecture.

0 1 0
̇ =[
x(𝑡) 𝑘 𝑑 ] x(𝑡) + [ 1 ] 𝑢(𝑡)
−𝑚 −𝑚 𝑚

While it is already stable we can still improve its performance by enhancing the dampening and
disturbance rejection. The controller should be a linear controller that takes into account the cur-

11.4. Optimal Control 145


Python for Engineers

rent deviation in position 𝑦 = 𝑥1 and velocity 𝑣 = 𝑥2 :

𝑢(𝑡) = −k𝑇 x(𝑡) = −𝑘1 𝑥1 (𝑡) − 𝑘2 𝑥2 (𝑡)

import control as ct

k = m = 1
d = 0.1

oscillator = [Link](
[Link]([[0, 1], [-k / m, -d / m]]),
[Link]([[0], [1 / m]]),
[Link](2),
[Link]((2, 1)),
)

closed_loop = [Link]([0.5, 0.5])


ct.impulse_response([oscillator, closed_loop], T=20).plot(label=["open␣
↪loop", "closed loop"]);

To optimize the response, we need to define a loss function. Choosing a loss function is in some
sense similar to programming, we need to translate our intuitive idea of good or bad behavior in
to an abstract mathematical form that the optimizer can work with. It is often the most difficult part
i solving these kinds of problems.
We want to minimize the both the deviation in position 𝑦[𝑛] as well as velocity 𝑣[𝑛] while also
penalizing large control inputs 𝑢[𝑛] to avoid an overly aggressive controller. Whenever we have
multiple competing goals like these, we need to prioritize between them by giving them different
weights. This way the optimizer can take the relative importance of each goal into account to
arrive at an overall good solution.
In the absence of any better ideas for how to penalize a deviation, it is typically a good idea to
square it. This often yields the desired results while giving the optimizer a well-behaved function
to work with. With this we choose our loss as:
𝐿(𝑘1 , 𝑘2 ) = ∑ (𝐶𝑦 𝑦[𝑛]2 + 𝐶𝑣 𝑣[𝑛]2 + 𝐶𝑢 𝑢[𝑛]2 )
𝑛

11.4. Optimal Control 146


Python for Engineers

Where 𝐶𝑦 , 𝐶𝑣 , 𝐶𝑢 are the weights for the position, velocity and control input respectively.

C_y = 1
C_v = 1
C_u = 0.1

def loss_func(params: [Link]) -> float:


k1, k2 = params
closed_loop = [Link]([k1, k2])
impulse_response = closed_loop.impulse_response(T=10, timepts_
↪num=1000)
positions = impulse_response.x[0]
velocities = impulse_response.x[1]
control_inputs = -k1 * positions - k2 * velocities
return C_y * [Link](positions**2) + C_v * [Link](velocities**2) + C_u␣
↪* [Link](control_inputs**2)

res = [Link](loss_func, x0=[0.5, 0.5])

closed_loop_opt = [Link](res.x)
ct.impulse_response([oscillator, closed_loop, closed_loop_opt], T=20).
↪plot(
label=["open loop", "closed loop", "optimized closed loop"]
);

This kind of optimization is actually a very common problem in control called a linear quadratic
regulator (LQR). The name derives from the fact that the system is linear, and the cost is quadratic.
There are special methods available like control’s lqr that solve these kinds of problems very
efficiently.

res_lqr, _, _ = [Link](oscillator, [Link]([C_y, C_v]), C_u)


print(res.x)
print(res_lqr[0])

[2.31736438 3.5933049 ]
[2.31662479 3.72664992]

11.4. Optimal Control 147


Python for Engineers

The results are not quite equal which is mostly due to numerical errors and the fact that we only
considered a finite time horizon as opposed to lqr which assumes an infinite time horizon.

11.4.2 Online Optimization

We can also go a step further. Instead of optimizing the controller ahead of time, we can also
optimize the individual control inputs online. For this let’s assume we have a desired output tra-
jectory 𝑦∗ (𝑡) that we want the system to follow. We want to find the controls 𝑢(𝑡) that minimize the
deviation from this trajectory:
𝑇
∫ (𝑦(𝑡) − 𝑦∗ (𝑡))2 𝑑𝑡
0

To make it a bit more interesting we will also want to limit the control range such that

−𝑢max ≤ 𝑢(𝑡) ≤ 𝑢max

and we only have a certain control budget


𝑇
∫ |𝑢(𝑡)|𝑑𝑡 ≤ 𝑈max
0

that we must not exceed.


To conform to these criteria we put bounds on the parameters and formulate a NonlinearCon-
straint for the control budget.

ts = [Link](0, 5, 50)
desired_response = [Link](ts - 2, 0, 1)
u_max = 1.5
budget = 3

[Link](ts, desired_response, label="desired response")

def loss_func(controls: [Link]) -> float:


response = oscillator.forced_response(ts, controls).x[0]
return [Link]((response - desired_response) ** 2)

constraint = [Link](lambda controls: np.


↪trapezoid([Link](controls), ts), lb=-[Link], ub=budget)
res = [Link](
loss_func,
x0=np.zeros_like(ts),
bounds=[(-u_max, u_max)] * len(ts),
constraints=[constraint],
options={"maxiter": 1000},
)
print(res)
controls = res.x
print("total action:", [Link]([Link](controls), ts))
response = oscillator.forced_response(ts, controls).x[0]

[Link](ts, response, label="optimized response", ls="--")


[Link](ts, controls, label="optimized controls")
[Link]();

11.4. Optimal Control 148


Python for Engineers

message: Optimization terminated successfully


success: True
status: 0
fun: 0.013628313950195184
x: [ 5.408e-04 -6.922e-05 ... -4.258e-04 6.636e-04]
nit: 107
jac: [-1.084e-04 -6.427e-05 ... -5.753e-04 -9.387e-05]
nfev: 5488
njev: 107
multipliers: [ 1.660e-02]
total action: 2.999999999953391

We can see that the resulting controls compensate for and exploit the resonance by anticipating
the initial acceleration and compensating for the ringing after reaching the setpoint. Overall, the
optimization result is very close to the desired trajectory while conforming to all our restrictions.
This way of optimizing the control inputs online is a very simple form of model predictive con-
trol (MPC). If you are doing something with MPC you should probably use a library like do-MPC
instead of implementing this by hand.

11.5 Summary

In this lecture we learned about:


• How to solve equations and systems of equations numerically.
• How to find minima and maxima of functions potentially subject to constraints.
• The caveats associated with local root finding and optimization.
• Curve fitting as a special case of optimization.
• Formulating (optimal control) objective functions

11.5. Summary 149


CHAPTER

TWELVE

PYTHON PROJECTS

Open in JupyterLite: Solved Blank


In all previous lectures and exercises, we have been exclusively working with Jupyter notebooks.
While they are an invaluable tool for simple demonstrations and “exploratory” programming, they
are more of a tool than an end product, and as such, they are not suitable for actual projects.
Proper projects are much more than a Jupyter notebook, they are spread across multiple files,
organized into directories, accompanied by documentation and other resources.
Python does not have a definite authority for build systems, leading people to organize their
projects in many different ways. In this chapter, we are going to focus on the basic concepts
and language features required for all of these approaches.
Note: Jupyter / IPython allows executing shell commands in notebooks by prefixing them with !.
So cells like:

!python3 my_script.py

are standins for something that you would execute in your shell (bash, PowerShell, …) like this:

> python3 my_script.py

12.1 Importing Code from Other Files

An important part of clean project organization is splitting the code into multiple source files called
modules. To make these modules work together, we need to be able to access functions and
other definitions from one module in another.

12.1.1 Files in the Same Directory

In the simplest case, we have two files / modules in the same directory, let’s call them [Link]
and [Link]. We want to use definitions from [Link] in [Link].

%%writefile my_module.py

class MyClass:
pass

def my_function():
(continues on next page)

150
Python for Engineers

(continued from previous page)


print("Hello from my_function!")

MY_CONSTANT = 42

print("this is executed on (first) import")

There are multiple ways to access these definitions in main:


The simplest way is to simply import the whole module. This will run the file and make all resulting
definitions available in a namespace with the module’s name:

%%writefile main_a.py

import my_module

obj = my_module.MyClass()
my_module.my_function()
print(my_module.MY_CONSTANT)

this is executed on (first) import


Hello from my_function!
42

We can also import it under an alias, just like we do with numpy and other common libraries:

%%writefile main_b.py

import my_module as mod

mod.my_function()

!python3 main_b.py

this is executed on (first) import


Hello from my_function!

Importing something with an alias can be confusing to the reader since it is no longer obvious
where the definitions come from, so only do it if it significantly improves readability or the alias is
a common one like with np for numpy.
It is also possible to only import parts of a module using a from ... import ... statement:

%%writefile main_c.py

from my_module import my_function


from my_module import MyClass, MY_CONSTANT # import multiple comma␣
↪separated definitions

my_function()

!python3 main_c.py

this is executed on (first) import


Hello from my_function!

12.1. Importing Code from Other Files 151


Python for Engineers

When we make multiple import calls to the same module, it is only executed once, making the
second import much faster, if the initial import requires expensive operations.
We can also import with an alias when using from imports. This can be useful for avoiding names-
pace conflicts when two modules define the same name. However, in these cases it is often better
to just use the fully qualified name like [Link]() to distinguish them.
%%writefile main_d.py

from my_module import my_function as func

func()

Finally, it is also possible to import everything from a module using the wildcard * character. This
should mostly be avoided however, since it is not immediately clear which definitions come from
which module if there are multiple wildcard imports.
%%writefile main_e.py

from my_module import *

my_function() # it is not obvious where this function is defined

12.1.2 Files in Different Directories

Things get somewhat more complicated when the files are located in different directories. Con-
sider the following directory structure:
src/
├── main/
│ ├── dir_1/
│ │ └── module_1.py
│ │ > def func_1():
│ │ > print("hello from func_1")
│ └── [Link]
│ > ...
├── module_2.py
│ > def func_2():
│ > print("hello from func_2")
└── dir_3/
└── module_3.py
> def func_3():
> print("hello from func_3")

There are three notable directions of imports here:


Importing func_1 into [Link]
This on is the simplest. If the file we want to import from is in a directory that is a sibling to the
current file, the directory is considered a module and the file in it a submodule so we can address
it with dot notation:
%%writefile src/main/main_a.py

from dir_1 import module_1


(continues on next page)

12.1. Importing Code from Other Files 152


Python for Engineers

(continued from previous page)


module_1.func_1()

from dir_1.module_1 import func_1 # import only the function


func_1()

!python3 src/main/main_a.py

hello from func_1


hello from func_1

We need to use the from syntax here since dir_1 is not a Python file and does not contain an
__init__.py file, more on that here. In general subdirectories are addressed as submodules
using a ., in most cases, this . behaves very similar to the path separator / in the file system.
Importing func_2 into [Link]
This one is not as simple, we need to address a parent of the current file. This is something that
Python does not normally allow. If we attempt it, we get an error message:

%%writefile src/main/main_b.py

## .. is the parent directory to this file


from .. import module_2

module_2.func_2()

!python3 src/main/main_b.py

Traceback (most recent call last):


File "/builds/python-for-engineers/course/build/book/chapters/12_
↪python_projects/src/main/main_b.py", line 3, in <module>
from .. import module_2
ImportError: attempted relative import with no known parent package

The interpreter seems to know what we want to do, but it does not allow us to do it. This is because
it always considers the file executed as the root of the project, and it does not allow accessing files
higher up file system hierarchy, as they are not considered part of the package.
The clean solution to this is to execute the file with a different root that contains the file we want
to import. This can be done by telling python to consider the directory as a module using the -m
flag. Since the file we want to execute is now considered a submodule of the package we replace
the path separator / with a dot . and leave out the .py extension:

!python3 -m [Link].main_b

hello from func_2

In some cases you may not want to do this, and you simply want to run the file as a script. One
“hack” to work around this, is to add the (parent) directory we want to import from to the list of
locations that Python searches for modules. Normally this contains the current file’s directory and
other directories that contain installed libraries.

12.1. Importing Code from Other Files 153


Python for Engineers

import sys

[Link]

['/home/vscode/.pyenv/versions/3.12.9/lib/[Link]',
'/home/vscode/.pyenv/versions/3.12.9/lib/python3.12',
'/home/vscode/.pyenv/versions/3.12.9/lib/python3.12/lib-dynload',
'',
'/builds/python-for-engineers/course/.venv/lib/python3.12/site-packages
↪']

After appending the parent src directory to the [Link] list, we import the module as if we
were importing from that parent path:

%%writefile src/main/main_b_hack.py

import sys

[Link]("src")

import module_2

module_2.func_2()

!python3 src/main/main_b_hack.py

hello from func_2

Note that the path we have to append depends on the working directory of the shell we are exe-
cuting the script from. This is not necessarily a bad thing since it is not usually a good idea to just
blindly add the parent directory of the current file to the search path. It is possible however, by
using some path manipulation tricks:

import sys
from pathlib import Path
[Link](str(Path(__file__).absolute().parent))
# file is the current files path, it is not defined in an interactive␣
↪session like Jupyter notebooks

Importing func_3 into [Link]


Importing from a file in a sibling directory is very similar to the previous case and also requires
running the file as a module or using path hacks. The only difference is that we now need to
address the file in the directory as a submodule like we did in the first case:

%%writefile src/main/main_c.py

from ..dir_3 import module_3

module_3.func_3()

!python3 -m [Link].main_c

12.1. Importing Code from Other Files 154


Python for Engineers

hello from func_3

Or with the path hack:

%%writefile src/main/main_c_hack.py

import sys

[Link]("src")

from dir_3 import module_3

module_3.func_3()

!python3 src/main/main_c_hack.py

hello from func_3

This should cover most cases that you will encounter in practice. For more details on importing
you can take a look at the official Python tutorial.

12.2 Making Scripts Importable

As discussed in the previous section, importing a file will execute all top level statements in it. So
if we want to import a function from a script, it will be executed potentially leading to unwanted
side effects or errors. To prevent this, we can include a check that determines whether the file is
being executed as a script and only execute the top level code if it is.

%%writefile my_script.py

def useful_function():
return 42

def main():
print("this is only printed when executed as a script", useful_
↪function())

# `__name__` is either "__main__" when executed as a script or the file␣


↪name when imported
print(f"this always printed ({__name__ = })")

if __name__ == "__main__":
# only execute this when this is the top level module, ie. it is run␣
↪as a script
main()

!python my_script.py

this always printed (__name__ = '__main__')


this is only printed when executed as a script 42

The __name__ variable is set

12.2. Making Scripts Importable 155


Python for Engineers

%%writefile other_script.py

import my_script

print("the answer is", my_script.useful_function())

!python3 other_script.py

this always printed (__name__ = 'my_script')


the answer is 42

The main use case for this is writing scripts or applications that can also be imported into other files
to allow for easier interoperability without having to resort to the (limited) command line interface
of the script.

12.3 Minimal Project Setup

So, given this freedom of organization, how should one structure a Python project? As mentioned
in the intro there is not universally agreed upon structure, but this minimal setup should be a
reasonable choice most simpler projects that are not intended to be distributed.
project_directory/
├── docs/
│ └── ... # project documentation
├── tests/
│ └── ... # testing code and resources
├── src/
│ └── package_name/
│ └── ... # source files and directories for submodules
├── [Link] # list of dependencies or [Link] etc. for␣
↪more advanced setup
└── [Link] # project description and setup instructions

The structure should be mostly self-explanatory:


The docs directory contains project documentation, this can either be a few simple markdown
files or a more complex setup that automatically generates documentation from docstrings and
other sources. The most common tool for this is Sphinx.
The tests directory contains tests that should verify that the code in src works as expected,
more details on this in later.
The actual source code is located in a directory with the module name inside the src directory.
It is also common to leave out this double nested structure and only have a src directory or a
directory with the package name in directly in the project root.
Most projects or packages depend on other packages, these are either specified in a
[Link] file or as part of other project configuration files in a [Link]
file. More on this later.
Finally, the [Link] file should contain a short description of the project and how to install /
use it. Repository hosting services like GitHub or GitLab will display this file in the project root,
so it should be structured like a landing page for the project with the most important information
immediately visible.

12.3. Minimal Project Setup 156


Python for Engineers

12.4 Managing Dependencies

There are many ways to manage dependencies for Python projects. The simplest one is to simply
list all required packages in a [Link] file places in the project root. The file contains
one package and optional version constraint per line:

numpy # any version of numpy, usually resolved to the latest version


matplotlib==3.10.3 # exactly version 3.10.3
scipy>=1.15.0 # at least version 1.15.0

These files are usually created semi-automatically. You first install the packages in your virtual
environment using pip install and then export the list of installed packages and their depen-
dencies to a file using the pip freeze command:

pip freeze > [Link]

Or with conda list if you are using conda:

conda list --export > [Link]

The environment can then be recreated somewhere else by reinstalling the same version of all
the packages listed in the requirements file:

python3 -m venv .venv # create a virtual environment


source .venv/bin/activate # activate the environment (linux/mac)
# ./.venv/Scripts/Activate.ps1 # activate the environment (windows)
pip install -r [Link]

Or the equivalent with conda:

conda create --name env-name --file [Link]

There are also more advanced and modern ways to manage dependencies like Poetry or Pixi.
These tools also take care of other project configuration like build systems, packaging, and pub-
lishing.

12.5 Testing

Testing is an important part of any robust project. These automated tests allow us to verify that
the code works as expected and that changes do not break existing functionality. A very conve-
nient framework for this is pytest. It provides a simple way to run tests and display helpful error
messages when tests fail.
Say we have the following module in our source directory:

%%writefile src/calendar_utils.py

def is_leap_year(year: int) -> bool:


return year % 4 == 0 # (incorrect)
# return year % 4 == 0 and (year % 100 != 0 or year % 400 == 0)

To test it, we create a corresponding file in the tests directory with the same name as the module,
but with a test_ prefix:

12.4. Managing Dependencies 157


Python for Engineers

%%writefile tests/test_calendar_utils.py

import sys
[Link]("src") # for sibling directory import
# can (and should) alternatively be defined in `pythonpath` in [Link]␣
↪or [Link]

import calendar_utils

def test_is_leap_year():
# use assertions with expected results
assert calendar_utils.is_leap_year(2020) == True
assert calendar_utils.is_leap_year(2021) == False
assert calendar_utils.is_leap_year(2100) == False
assert calendar_utils.is_leap_year(2000) == True

To run the tests, we simply invoke pytest in the project root. It will then automatically discover
all tests in the tests directory and run them. It is usually a good idea to add the -v flag to make
the output more verbose providing a detailed report of the test results.

!pytest -v

============================= test session starts␣


↪==============================
platform linux -- Python 3.12.9, pytest-8.4.1, pluggy-1.6.0 -- /builds/
↪python-for-engineers/course/.venv/bin/python
cachedir: .pytest_cache
rootdir: /builds/python-for-engineers/course
configfile: [Link]
plugins: anyio-4.9.0
collecting ...

collected 1 item ␣

tests/test_calendar_utils.py::test_is_leap_year

FAILED [100%]

=================================== FAILURES␣
↪===================================
______________________________ test_is_leap_year _______________________
↪________

def test_is_leap_year():
# use assertions with expected results
assert calendar_utils.is_leap_year(2020) == True
assert calendar_utils.is_leap_year(2021) == False
> assert calendar_utils.is_leap_year(2100) == False
E assert True == False
E + where True = <function is_leap_year at 0x78969af6f6a0>(2100)
E + where <function is_leap_year at 0x78969af6f6a0> =␣
↪calendar_utils.is_leap_year

12.5. Testing 158


Python for Engineers

tests/test_calendar_utils.py:12: AssertionError
=========================== short test summary info␣
↪============================
FAILED tests/test_calendar_utils.py::test_is_leap_year - assert True ==␣
↪False
+ where True = <function is_leap_year at 0x78969af6f6a0>(2100)
+ where <function is_leap_year at 0x78969af6f6a0> = calendar_utils.
↪is_leap_year

============================== 1 failed in 0.12s␣


↪===============================

The test fails because the is_leap_year() function is not implemented properly. The test
report shows us exactly which assertion failed and how it differs from the expected value. Years
that are divisible by 100 are not leap years unless they are also divisible by 400. After adding
these additional rules, all tests should pass successfully.

%%writefile src/calendar_utils.py

def is_leap_year(year: int) -> bool:


return year % 4 == 0 and (year % 100 != 0 or year % 400 == 0)

!pytest -v

============================= test session starts␣


↪==============================
platform linux -- Python 3.12.9, pytest-8.4.1, pluggy-1.6.0 -- /builds/
↪python-for-engineers/course/.venv/bin/python
cachedir: .pytest_cache
rootdir: /builds/python-for-engineers/course
configfile: [Link]
plugins: anyio-4.9.0
collecting ...
collected 1 item ␣

tests/test_calendar_utils.py::test_is_leap_year PASSED ␣
↪ [100%]

============================== 1 passed in 0.03s␣


↪===============================

Writing good tests is often not easy and sometimes only possible in hindsight. Ideally the tests
should cover every possible path through the control flow, but this is rarely feasible. A (somewhat)
more realistic goal is to make sure that the tests cause every line of code to be executed at least
once. This is called code coverage and can be measured using tools like [Link].
Testing should not be overdone, however. There is no point in writing 10 tests that all test the
same code path just with slightly different inputs. Write a test for every edge case or special case
(if feasible), and a few tests for the common cases. If you fixed or at least identified a bug, write a
test that fails because of that bug to avoid reintroducing it in the future. Focus on writing tests for
features, that are hard to get right on the first try and end-to-end tests that verify single use cases
from start to finish.

12.5. Testing 159


Python for Engineers

12.6 Summary

In this lecture we learned about:


• How to import code from other files and directories.
• Basic Python project structure.
• Managing dependencies with a [Link] file.
• Writing and running tests with pytest.
Some of the things we did not cover are:
• Version control with git and hosting services like GitHub or GitLab.
• Writing and automatically generating documentation with Sphinx.
• Advanced dependency management, packaging, and distribution with Poetry or Pixi.
• And many other tools involved in developing and maintaining larger projects…

12.6. Summary 160

You might also like