Research Computing in Python

#1 Fundamentals — a short history of computing, and a why, what and how of Python

2026-07-17

Computing and Me

  • Dr Matthew Evans (https://ml-evs.science, https://github.com/ml-evs)
  • Research Fellow in the Grey Group, Department of Chemistry
  • Background in computational physics & materials science
  • Running hundreds of thousands of HPC jobs simulating materials with DFT
  • Now writes open source data managemement software (datalab) for materials and chemistry

Peak: #29 GitHub user in Belgium

This course

  • A biased take on what you need to get started using computers in research
  • Not a computer science course, not exhaustive
  • A mixture of high-level concepts, empirical truths and practical advice
  • Not a real Python course, either - you can’t learn a language by reading the dictionary

This lecture: Fundamentals

  • Some computing history
  • Some computing architecture
  • Programming languages
  • Data structures
  • Python
    • Why Python is slow — and why we use it anyway
    • A demo: the Mandelbrot set
  • Wrap-up and motivation for next week

A brief history of computing

A brief history of computing

The Antikythera mechanism

Babbage & Lovelace engine

Images: Wikimedia Commons — Antikythera (CC BY-SA 3.0), Difference Engine (public domain)

Douglas Hartree and his differential analyser

Meccano differential analyser for computing atomic spectra

\(1\,E_\mathrm{h} = \dfrac{\hbar^2}{m_\mathrm{e} a_0^2} = \dfrac{e^2}{4\pi\varepsilon_0 a_0} \approx 27.2\ \mathrm{eV}\)

Images: differential analyser (Sylvain Machefert, CC BY-SA 4.0); portrait (fair use)

“The implications of the machine are so vast that we cannot conceive how they will affect our civilisation. Here you have something which is making one field of human activity 1,000 times faster. In the field of transportation, the equivalent to ACE would be the ability to travel from London to Cambridge … in five seconds as a regular thing. It is almost unimaginable.”

— Douglas Hartree, on the Automatic Computing Engine, 1946

…and yet, all the calculations that would ever be needed in this country could be done on three digital computers: one in Cambridge, one in Teddington, and one in Manchester.

— Hartree in 1951, as recalled by Lord Bowden

The Manchester “Baby” (1948)

EDSAC, Cambridge (1949)

Images: Manchester Baby replica (Logg Tandy, CC BY 4.0); EDSAC (University of Cambridge Computer Laboratory, CC BY 2.0)

The first stored-program computers - the birth of software?

The first actual computer bug (Harvard Mark II, 1947)

Image: Wikimedia Commons (public domain, US Navy)

Parts of the first four US Army computers (ENIAC, EDVAC, ORDVAC, BRLESC), 1945–1962

BRLESC-I in full

Images: Wikimedia Commons — Army computers, BRLESC-I (both public domain, US Army)

A brief history of computing

The first transistor (1947)

1 transistor

Intel 4004 (1971)

2,300 transistors

Apple M1 (2020)

16,000,000,000 transistors

In the meantime… hard disks, programming languages, operating systems, parallel computers, graphics, personal computers, networking, the Internet, the World Wide Web

Images: Wikimedia Commons — first transistor (public domain), Intel 4004 (Thomas Nguyen, CC BY-SA 4.0), Apple M1 (CC0)

What have computers been used for?

  • First, science — tables, trajectories, atomic structure: machines built by and for researchers
  • Then, the military: ballistics, cryptography, bomb guidance
  • Then, commerce: payroll, banking, inventory
  • Then, suddenly, everyone: PC revolution + the internet/web democratised computing: communication, creativity, everyday life
  • Yet the frontier is still pushed by science and nowadays by AI

LUMI, Kajaani, Finland (2022)

MareNostrum 4, Barcelona (2017)

Images: LUMI (LUMI consortium / Fade Creative); MareNostrum 4 (Gemmaribasmaspoch, CC BY-SA 4.0)

Software in research today

  • Simulations: molecular dynamics, quantum chemistry and materials science, climate, weather, astrophysics, complex systems
  • Automation: instrument control, robotics, high-throughput experiments, laboratory information management systems
  • Data analysis: large-scale experiments, imaging, plotting
  • Statistical & probabilistic modelling: fitting and inference, uncertainty quantification, Monte Carlo and Bayesian methods
  • Data management: provenance, reproducibility, FAIR principles
  • Machine learning: image recognition, natural language processing, generative models
  • Dissemination: publishing, visualisation, outreach, citizen science

Computing architecture

What is the machine actually doing?

Intel 8008 (1972), die labelled

  • The CPU receives instructions as “machine code”, human-readable as assembly
  • decodes each one
  • and executes it, shuttling values between the registers (right) and the arithmetic/logic unit (left)
  • All chip design since revolves around efficiently cramming as much compute as possible into this layout: specialised registers, specialised instructions…
  • …and keeping it fed: modern CPUs pipeline instructions, perform branch prediction execute speculatively and out of order, and hide slow memory behind layers of caches

Image: Ken Shirriff, righto.com — used with attribution

Instruction sets: x86 and ARM

  • x86: the Intel 8086 (1978) instruction set — descendant of the 8008, still running nearly every “Intel or AMD” machine
  • amd64 / x86-64: its 64-bit extension — designed by AMD (2003), adopted by Intel, extended to include vectorisation to make the most of wider registers (SSE, AVX)
  • ARM (Acorn Computers, Cambridge, 1985): the only real competitor — every phone, Apple M-series, now Windows laptops

Programming

Programming languages

Programming languages are the “simple”, human-facing interface to all this machinery.

What you write (C):

double mean(const double *x, int n) {
    double total = 0.0;
    for (int i = 0; i < n; i++)
        total += x[i];
    return total / n;
}

What the compiler turns it into (x86-64, gcc -O1):

mean:
        test    esi, esi
        jle     .L4
        mov     rax, rdi
        movsxd  rdx, esi
        lea     rdx, [rdi+rdx*8]
        pxor    xmm0, xmm0
.L3:
        addsd   xmm0, QWORD PTR [rax]
        add     rax, 8
        cmp     rax, rdx
        jne     .L3
.L2:
        pxor    xmm1, xmm1
        cvtsi2sd xmm1, esi
        divsd   xmm0, xmm1
        ret
.L4:
        pxor    xmm0, xmm0
        jmp     .L2

Low-level vs high-level languages

  • The “level” of a language ≈ its distance from the machine
  • C (1972): low-level — you manage the memory, and the compiler maps your code almost 1:1 onto the hardware
  • Fortran (1957): the original high-level language — “FORmula TRANslation”, built for science; compiled and fast, with arrays as first-class citizens
  • Python (1991): very high level — no types to declare, memory managed for you, everything decided at runtime… which is exactly why it’s slow
/* the max of some numbers, in C */
double *xs = malloc(5 * sizeof(double));
if (xs == NULL)                /* did we even get the memory? */
    exit(1);
xs[0] = 1.0; xs[1] = 4.5; xs[2] = 2.0;
xs[3] = 8.1; xs[4] = 3.0;

double largest = xs[0];
for (int i = 1; i < 5; i++)
    if (xs[i] > largest)
        largest = xs[i];

free(xs);                      /* give it back, or leak it */
! the max of some numbers, in Fortran
real :: xs(5) = [1.0, 4.5, 2.0, 8.1, 3.0]
real :: largest
largest = maxval(xs)   ! arrays built in
# the max of some numbers, in Python
largest = max([1.0, 4.5, 2.0, 8.1, 3.0])

Software usage on ARCHER2, the UK national supercomputer (June 2026)

Data: ARCHER2 usage data, EPCC

A programming language is a contract

  • The language is abstract: syntax and rules of meaning, written down in a standard, e.g., the C and Fortran standards, the Python language reference. No machine is mentioned
  • The implementation honours that contract on a particular machine: gcc/clang for C, gfortran/flang for Fortran, CPython/PyPy for Python
  • The mean() you just saw compiles to different assembly on ARM but the C doesn’t change. That’s what portable means: write to the contract, not to the machine

Everything is an abstraction

        your analysis script
              Python
      NumPy / compiled libraries
          compiler ・ OS
        assembly ・ x86, ARM
        CPU ・ gates ・ transistors
              electrons
  • Computing is a tower of abstractions — nobody holds it all in their head at once
  • The working rule: know enough to reason one layer down, and to design for the layer above
  • But abstractions leak — the Law of Leaky Abstractions: the layer below pokes through exactly when things go wrong
  • 0.1 + 0.2 ≠ 0.3, “why is this loop slow?”, “why won’t this install on my machine?”

Data structures

Data structures: bits and integers

Everything in memory is bits — a type is a rule for reading them.

The simplest type: a boolean, one bit of information, true or false:

bool done = false;
0

A 32-bit signed integer:

int n = 42;    /* exactly 4 bytes */
00000000 00000000 00000000 00101010

Range: \(-2^{31}\) to \(2^{31}-1\) = 2,147,483,647 — one more and it overflows:

int n = 2147483647 + 1;
printf("n = %d\n", n);
n = -2147483648
unsigned int u = 2147483648 + 1;
printf("u = %u\n", u);
u = 2147483649

The same 32 bits, read as a float:

float x = 42.0f;   /* also 4 bytes */
01000010 00101000 00000000 00000000
[ sign | 8-bit exponent | 23-bit fraction ]

Range: \(\pm 3.4 \times 10^{38}\) (smallest positive \(\approx 1.2 \times 10^{-38}\)) — but only ~7 significant digits

  • Y2038: signed 32-bit Unix time (number of seconds since 01/01/1970) runs out at 03:14:07 UTC on 19 January 2038

Data structures: floating point (double precision)

A double is 64 bits read as \(\pm\, m \times 2^{e}\) (IEEE 754); used extensively in scientific computing.

[ 1 sign bit | 11 exponent bits | 52 fraction bits ]
double x = 0.1;               /* stored: 0.10000000000000000555… */
printf("%.17f\n", 0.1 + 0.2);   /* prints 0.30000000000000004 */
  • 0.1 has no exact binary representation — like 1/3 in decimal
  • so never compare floats with ==; test fabs(a - b) < tol
  • smaller floats (32-, 16-bit) trade precision for speed — the currency of GPUs and machine learning

Arrays and vectors

double xs[1024];              /* C: one contiguous block, size fixed forever */
std::vector<double> ys = {1.0, 4.5};  // C++: contiguous but growable
ys.push_back(2.0);   // may silently reallocate-and-copy the whole thing
  • Contiguity is why arrays are fast: caches and vector units love them

Python

The locker and the wardrobe

A C variable is a locker:

int x = 42;      /* 4 bytes at a fixed address */
x = "hello";     /* compile error: wrong shape */
  • fixed size, fixed shape, the name is the address

A Python object is the wardrobe to Narnia:

x = 42          # a name, tied to an object
x = "hello"     # fine — retie the name
  • looks like a small box, but it’s bigger on the inside: the object carries its type, its size, its reference count…

Image: another Nano Banana hallucination

Why “Python is slow”

  • Every value is an object on the heap — even 1 (a C int is 4 bytes; a Python int is 28)
  • Everything is mutable, so the interpreter can assume nothing, so cannot optimise anything ahead of time
  • Types are checked at runtime, on every operation, every time
  • A list is an array of pointers to scattered objects, no cache locality

…so why does everyone use it?

  • Clean syntax: executable pseudocode, optimised for reading
  • Hides complexity: memory management, types, compilation, linking, installation
  • The ecosystem: “batteries included” standard library, plus a mature external library for almost anything (PyPI)
  • Offloading: NumPy, PyTorch & friends are compiled C/C++/Fortran underneath
  • Economics: researcher time usually costs more than CPU time

Higher-level structures: the dictionary

In C, you build it yourself:

struct entry {
    const char *key;
    double value;
    struct entry *next;
};
unsigned hash(const char *s);   /* you write this */
double lookup(struct entry **t,
              const char *key); /* and this */
/* …and collisions, resizing, deletion,
   and who owns all that memory? */

In Fortran: no dictionary in the standard library at all.

In Python, it’s a literal:

masses = {"H": 1.008, "He": 4.003}
masses["Li"] = 6.94

Python is built on dicts: module namespaces, object attributes, keyword arguments — dicts all the way down.

Dicts speak JSON

A Python dict…

sample = {
    "formula": "NaCl",
    "mass_g": 0.253,
    "conductivity": [1.2, 1.4, 1.1],
}

import json
json.dump(sample, open("sample.json", "w"))

…is one line away from a file anything can read:

{
  "formula": "NaCl",
  "mass_g": 0.253,
  "conductivity": [1.2, 1.4, 1.1]
}
  • JSON is just dicts, lists, numbers, strings and booleans — the lingua franca of data exchange
  • Generic formats outlive any one program, language, or PhD: this file is readable from Python, R, Julia, C++, a browser, a database — today and in twenty years
  • Your instrument’s proprietary binary format is none of those things
  • Rule of thumb: keep open, generic formats (JSON, CSV, HDF5, …) at the boundaries of everything you write

The Python types you’ll actually use

Type Literal Notes
bool True, False
int 42 arbitrary precision — no overflow!
float 6.626e-34 an IEEE double underneath
str "NaCl" Unicode text
list [1, 2, 3] mutable sequence
tuple (x, y) immutable sequence
dict {"H": 1.008} key → value
set {2, 3, 5} uniqueness & membership
None None “nothing to see here”
class class Atom: ... user-defined objects

Example: Periodic Table, four ways

  1. every element a variable?
H = 1
He = 2
Li = 3           # …116 more lines of this
print(H)                     # 1
  1. How about a list?
elements = ["H", "He", "Li", "Be", "B", "C", "N", "O", ...]
elements.index("C") + 1      # 6 — mind the off-by-one!
elements[8+1]                    # "O"

Loopable, and lookup works — but one fact per element, hidden in its position

  1. a dictionary: keyed by what you actually know
ptable = {
    "H":  {"Z": 1, "mass": 1.008, "name": "hydrogen"},
    "He": {"Z": 2, "mass": 4.003, "name": "helium"},
}
ptable["He"]["mass"]       # 4.003

Look things up by meaning, attach as many properties as you like

Example: Periodic Table, four ways

  1. A custom class
class Element:
    def __init__(self, symbol, Z, mass, name):
        self.symbol = symbol
        self.Z = Z
        self.mass = mass
        self.name = name

ptable = {
    "H":  Element("H", 1, 1.008, "hydrogen"),
    "He": Element("He", 2, 4.003, "helium"),
}

Control flow in Python

Loops (for, while), conditionals (if, elif, else) and functions build programs

ptable = {
    "H": Element("H", 1, 1.008, "hydrogen"),
    "He": Element("He", 2, 4.003, "helium"),
    ...
}

def is_element(symbol):
    """Is the symbol a known element?"""
    if symbol in ptable:
        return True
    return False

is_element("H")  # True
is_element("Xx") # False
is_element(8) # False

Control flow in Python

Loops (for, while), conditionals (if, elif, else) and functions build programs

ptable = {
    "H": Element("H", 1, 1.008, "hydrogen"),
    "He": Element("He", 2, 4.003, "helium"),
    ...
}

def is_element(symbol):
    """Is the symbol a known element?"""
    return symbol in ptable

is_element("H")  # True
is_element("Xx") # False
is_element(8) # False

Control flow in Python

Loops (for, while), conditionals (if, elif, else) and functions build programs

ptable = {
    "H": Element("H", 1, 1.008, "hydrogen"),
    "He": Element("He", 2, 4.003, "helium"),
    ...
}

def is_lanthanide(symbol):
    """Is the element a lanthanide?"""
    if symbol in ptable and 57 <= ptable[symbol].Z <= 71:
        return True

Control flow in Python

Loops (for, while), conditionals (if, elif, else) and functions build programs

ptable = {
    "H": Element("H", 1, 1.008, "hydrogen"),
    "He": Element("He", 2, 4.003, "helium"),
    ...
}

def is_lanthanide(symbol):
    """Is the element a lanthanide?"""
    if symbol in ptable and 57 <= ptable[symbol].Z <= 71:
        return True

def get_lanthanides():
    """Return a list of all lanthanides"""
    lanthanides = []
    for symbol in ptable:
        if is_lanthanide(symbol):
            lanthanides.append(ptable[symbol])

    return lanthanides
  • Indentation is the syntax — blocks are visual
  • Express logical combinations with operators and, or, not
  • Conditionals are expressions that can be chained (if, elif, else)
  • for iterates over any container (lists, dicts, files, …)

Control flow in Python

Loops (for, while), conditionals (if, elif, else) and functions build programs

ptable = {
    "H": Element("H", 1, 1.008, "hydrogen"),
    "He": Element("He", 2, 4.003, "helium"),
    ...
}

def is_lanthanide(symbol):
    """Is the element a lanthanide?"""
    if symbol in ptable and 57 <= ptable[symbol].Z <= 71:
        return True

def get_lanthanides():
    """Return a list of all lanthanides"""
    lanthanides = []
    for symbol in ptable:
        if is_lanthanide(symbol):
            lanthanides.append(ptable[symbol])

    return lanthanides

la_mass_sum = sum(e.mass for e in get_lanthanides())

Scientific Python

Nearly all of it stands on one shared abstraction: the NumPy array.

  • NumPy — the n-dimensional array: memory laid out like Fortran, driven from Python
  • SciPy — the numerical toolbox: integration, optimisation, linear algebra, statistics, signal processing
  • matplotlib — plotting (most of the figures in most of the papers you’ll read)
  • pandas — labelled tables (“dataframes”) for messy real-world data
  • Jupyter — notebooks: code, figures and narrative in one document
  • scikit-learn — classical machine learning; PyTorch / JAX at the deep-learning end
  • …and domain layers on top: astropy (astronomy), ASE / pymatgen (materials), RDKit (molecules), Biopython, MDAnalysis, …

NumPy: putting the speed back

Remember why arrays are fast — one contiguous block. NumPy brings that to Python:

import numpy as np

measurements = [[1.2, 1.4, 1.1],     # a nested list of Python objects…
                [0.9, 1.3, 1.5]]
arr = np.array(measurements)         # …becomes one contiguous block of float64
2 * arr             # array([[2.4, 2.8, 2.2],
                    #        [1.8, 2.6, 3. ]])
arr.mean()          # 1.2333333333333332
arr.max(axis=1)     # array([1.4, 1.5]) — the max of each row

One line of Python = one compiled loop over the whole block

arr[0]              # first row
arr[:, 2]           # third column: array([1.1, 1.5])
arr[arr > 1.2]      # by condition:  array([1.4, 1.3, 1.5])

Slicing and masking replace most of the loops you’d otherwise write

Optional final demo: Mandelbrot

The Mandelbrot set

Iterate one very simple map:

\[z_{n+1} = z_n^2 + c, \qquad z_0 = 0\]

for a complex number \(c\) — i.e. a point in the plane.

  • The set is every \(c\) for which \(z_n\) stays bounded forever
  • We can’t iterate forever, but there’s an escape test: once \(|z_n| > 2\), the orbit is gone for good — so iterate each pixel’s \(c\) up to max_iter times and colour it by when it escaped
  • The boundary is a fractal: infinitely detailed, new structure at every zoom
  • And it’s a perfect benchmark: millions of pixels × hundreds of iterations of pure arithmetic, nothing to hide behind

Demo: Mandelbrot

Pure Python:

x = y = 0.0
n = 0
while (x*x + y*y < 4.0
       and n < max_iter):
    x, y = (x*x - y*y + cx,
            2.0*x*y + cy)
    n += 1
out[j][i] = n

NumPy (whole grid at once):

for n in range(max_iter):
    z[mask] = z[mask]**2 + c[mask]
    escaped = mask & (np.abs(z) > 2)
    out[escaped] = n
    mask &= ~escaped

Fortran (via f2py):

do while (x*x + y*y < 4.0d0 &
          .and. n < max_iter)
    xt = x*x - y*y + cx
    y  = 2.0d0*x*y + cy
    x  = xt
    n  = n + 1
end do
Pure Python NumPy Numba Fortran
1024², 100 iterations 1.69 s 0.3 s 0.02 s 0.01 s

Wrap-up

Recap

  • A century of hardware, one enduring idea: the stored program, fetch–decode–execute — everything else is miniaturisation
  • Languages are abstractions over that machine: compiled ones (C, Fortran) stay close to it; Python trades machine time for your time
  • Types are rules for reading bits — and every abstraction leaks (overflow, 0.1 + 0.2, “why is this slow?”)
  • Python hides the machinery in objects and dicts; when you need speed, offload to the compiled layer (NumPy) — or fix the algorithm first

Not covered

  • Randomness: pseudo-random numbers, seeds, and why your “random” results should be reproducible
  • Recursion: functions that call themselves
  • Classes and object-oriented design — you’ll meet them the moment you use a library
  • Parallelism: threads, processes, MPI, GPUs
  • Exceptions and error handling
  • Algorithmic complexity, properly (big-O beyond today’s hand-waving)
  • Text encodings, regular expressions, databases, the network…

Software is a living ecosystem

  • Everything you’ve seen today changes: languages, libraries, hardware — software grows organically around historical biases (x86’s baggage, Fortran’s arrays, Python’s dicts)
  • You are part of the ecosystem: the code you write and the tools you choose feed back into it
  • So the same code will not always work, everywhere, forever — unless you fence it in: pinned dependencies, recorded environments, containerisation

Next time: programming in practice

  • Practical software development: git, GitHub, and collaborating on code
  • Research software engineering: best practices — and the anti-patterns they replace
  • Environments and packaging (uv), testing, documentation
  • Some thoughts on AI code generation and agents: what changes, what doesn’t

Further reading