Motivation ---------- I would like to start this tutorial for pyoomph with a little history to shed some light into the decision to develop pyoomph. Feel free to skip this section! In 2015, I started my postdoc at the TU Eindhoven with Professor Hans Kuerten where I generalized an existing Fortran code for the evaporation of sessile droplets to account for multi-component mixture droplets :cite:`diddens2017modeling,tan2016evaporation`. Since this code was based on lubrication theory, a reviewer of my first publication on this topic :cite:`diddens2017modeling` asked for a validation at higher contact angles, i.e. where the lubrication theory might fall short. Inspired by the work of Hu & Larson on pure droplets :cite:`Hu2005FEM,Hu2005Marangoni`, who used a finite element method for this validation step, I found FEniCS (`fenicsproject.org `_) :cite:`FEniCS1,FEniCS2` as a quick and simple solution to easily use the finite element method, but only after I noticed that writing a finite element library myself is quite demanding. It is in fact reinventing the wheel given all the nice existing open source libraries. With FEniCS, things are so simple: just type a few lines in Python and you directly get your solution - perfect! This nice experience gave rise to move the previous lubrication theory Fortran code to a full Navier-Stokes simulation of evaporating multi-component droplets in FEniCS :cite:`diddens2017detailed,tan2017self,tan2022evaporation,diddens2017evaporating,li2019gravitational,li2018evaporation`. However, while FEniCS is very powerful in quickly developing simple simulations, it has serious downsides for free surface problems and multi-physics scenarios, where individual equations have to be solved on different domains. Also the arbitrary Lagrangian-Eulerian method (ALE) cannot, to my best knowledge, be well implemented in FEniCS, at least not in a monolithic solution approach. This method, however, is very powerful for the problems I was working on that times. By accident, I came across the finite element library oomph-lib (`oomph-lib.maths.man.ac.uk `_) :cite:`Heil2006`. This toolbox offered everything I was looking for: free surfaces, multi-physics, ALE methods and the possibility to add surfactants were already ready to be used more or less out of the box. Even more, spatial and temporal adaptivity, arc length continuation and bifurcation tracking are already part of oomph-lib. So I again rewrote the code, now with the help of the very flexible and powerful library oomph-lib. In direct comparison with FEniCS, however, the definition of custom equations in oomph-lib is a demanding and time-consuming task, since the entire assembly of the spatially discretized weak forms of your equations have to be written by hand in C++. In fact, for best speed and convergence behavior, it is even required to also fill the Jacobian matrix of the weak form by hand, which easily ends up in writing 1000 lines of code for a rather simple equation. While my work with oomph-lib resulted in several publications :cite:`li2019bouncing,gauthier2019self,li2020evaporating,diddens2021competing,cleve2021time,hack2021asymmetric,de2021marangoni`, the limitations became more and more obvious: Each new equation you want to introduce in your system requires a lot of code writing and compilation, as well as substantial code overhead to couple all equations at the end. I was longing to return to the quick development approach I had with FEniCS, but getting all the benefits of oomph-lib. I therefore decided to aim for the best of both worlds - and pyoomph was born\ :math:`\ldots` In 2022, Duarte Rocha joined our group as PhD student and started to contribute to the development of pyoomph.