Preview

CITCOM is a geodynamics modelling code based on the finite element method that is designed for planetary evolution modelling where large spatial variations and strong non-linearities occur in the material properties.

Underworld3 has about 50,000 lines of Python/ Cython wrapping PETSc, SymPy, and a just-in-time compiler. I began a trial of AI coding tools in 2025 and they have gradually become central to the way our team works. This is a story of co-evolution as much as it is about adoption of a new set of tools.

A constitutive model is a Python class where the relationship between fluxes and gradients is encoded in SymPy. At every stage the mathematics is visible, inspectable, and differentiable. The framework handles Jacobians, C code generation, and PETSc integration. You handle the physics.

It can be difficult to figure out where you are in an unstructured mesh of triangles or tetrahedra; worse when the mesh is distributed in parallel. We have to solve this problem for particle-in-cell type codes such as Underworld3. This is how we do it.

As the old truism states, the best way to learn is by doing. This guide can get you from zero to hero in 60 seconds - give or take.

Parallel computation puts many CPUs to work on solving a problem much more quickly than one CPU alone. But this only works if the tasks are carefully scheduled and the additional CPUs are not waiting around for something to do. How do we choose the right number of processors for a given problem ?

The new generation of Apple Mac comes with the new Apple Silicon (M1) chip which has an Arm architecture (as opposed to the older generation that had i386 Intel processor). This brings all manner of troubles and requirements for the development of codes.

What actually happens between the moment you write a mathematical expression in underworld3's python layer and the moment PETSc receives a finite element term in the form of compiled C code ?

One of the less glamorous but important problems in a finite element framework is this: how does the user assign values to a field variable, and how does the framework ensure that PETSc sees those values correctly, in parallel, without the user needing to keep track themselves?

There are three ways to put resolution where a calculation needs it: rebuild the mesh, subdivide it, or move its nodes. The third conserves the node budget and leaves the connectivity, the parallel partition and the multigrid hierarchy untouched — nothing is created or destroyed, only relocated. How that redistribution works in Underworld3, what it is worth in numbers, and where its ceiling is.

Underworld is a code for geodynamics — mantle convection, lithospheric deformation, subduction, ice flow. We solve coupled, nonlinear PDEs with complex rheologies using Lagrangian particles to track material history. The project has been running for 20 years: why did we start again ?

Underworld is built around the idea of active Lagrangian tracer particles that carry history and composition information as the material deforms. How do we combine this information with our symbolic mathematical framework ?

Non-dimensionalisation (rewriting problems to make them re-scalable from lab to the real world) is a specialised task that all modellers confront at some point. It is important for accurate and efficient numerical solutions. Can we take away the pain that comes with the task ?

One link opens Underworld in a browser, with any public repository pulled in beside it and any released version underneath. Nothing is installed, and nothing has to be added to the repository being launched. We use it for teaching, for workshops, and for anyone who would rather run the model from a paper than read about it.

To test scalability we run weak scaling tests on various HPC machines to check the numerical framework remains robust when pushing for higher fidelity models.

How we built a simple dashboard using Github actions with open source software and openly available (FAIR) data.

Underworld can precondition a Stokes solve with geometric multigrid built from a real hierarchy of meshes, rather than with the algebraic multigrid it falls back to. What the difference is, how to build a mesh that has a hierarchy, and when the choice is worth making.

In Underworld3, the time derivative is a symbolic object. It appears in the solver's strong form as a SymPy expression, alongside the constitutive stress and the body force.

TL; DR: Underworld is now Singularity-enabled, making it easier and somewhat quicker to use compared to traditional HPC installs. We demonstrate results using over 10,000 CPUs and more than one billion unknowns for solving the Stokes equation with Underworld 2.16

We built a cheap-and-cheerful solution with persistent storage and a binder-like access to notebooks in a repository that is aimed at serving a single classroom. The zero-to-server time is just a few minutes and there is minimal manual configuration.

Introducing Underworld3: Mathematically Self-Describing Modelling in Python for Desktop, HPC and Cloud.