## 2 Mathematical Modelling and Numerical Methods
### Pressure Correction and Poisson Solvers  
In the second fractional step, the core is to first solve the pressure Poisson equation using an efficient numerical method and then apply the obtained pressure gradient to correct the intermediate velocity according to Equation (4) to enforce a divergence-free velocity field. Conventionally, the pressure Poisson equation is solved in the same computational domain as for the momentum equations, although the grid points for pressure can be different from those for velocity (e.g., a staggered grid in which the pressure variable is stored at the cell centres of the control volumes, whereas the velocity variables are located at the cell faces). There are many methods available for the numerical solution of the pressure Poisson equation. In this study, we adopted the WLS-based collocation method [31] as a reference method. Basically, in this method, the second-order derivatives in Equation (5) and the first-order derivatives in the boundary conditions are approximated at each collocation point by the weighted least-squares method described in the previous subsection. The Dirichlet boundary condition is satisfied by simply prescribing the pressures on the corresponding boundary points to the fixed values. If there is no Dirichlet boundary condition, the pressure at a chosen point should be fixed to a reference pressure value. In [31], the resulting sparse linear algebraic equations for the unknown pressure values at the present time level \(n\) were solved by an iterative method known as the preconditioned biconjugate gradient method [33]. In terms of time to solution, this is acceptable for small- to medium-sized linear systems, but too slow for large linear systems. Hence, it was replaced here by an algebraic multigrid method [34; 35; 36] implemented in the open-source program AGMG developed by [PERSON] (see [[http://agmg.eu](http://agmg.eu)]([http://agmg.eu](http://agmg.eu)) (accessed on 1 June 2023) for detailed documentation). The iteration starts with initial pressure values taken as those from the previous time level (note: the initial pressure field at time \(t=0\) should be given). As an alternative to conventional ways of performing the pressure-correction step, we propose a hybrid algorithm that solves the pressure Poisson equation in an extended and rectangular computational domain and performs the pressure correction to the velocity field in the flow domain as usual. The advantage of the hybrid algorithm is that it enables the use of the fast Fourier transform (FFT) for directly solving the Poisson equation on a uniform orthogonal grid, hence has the potential to speed up the overall computation in the second fractional step. However, it is important to note that it inherits the two main disadvantages of FFT-based methods: 1. The reliance on regular grids makes it less efficient for flow problems discretised with highly heterogeneous grid resolutions (e.g., due to local grid refinement). 2. It only works for certain pressure boundary conditions (to be discussed later). As illustrated in Figure 1, the pressure Poisson equation is solved in a rectangular domain (shown here as a square for convenience) discretised by uniformly distributed points (marked as red in the figure and, hereafter, called the pressure points) and the real flow domain is discretised by irregularly distributed points (marked as blue in the figure and, hereafter, called the velocity points). Like in the immersed boundary method [37], a solid obstacle (here, a cylinder marked by the black circle with its boundary discretised by the boundary points in black) is treated as an immersed body filled with a frozen fluid; hence, the flow field extends effectively to the whole rectangular domain with the incompressibility holding everywhere and the no-slip condition enforced at the immersed boundary (IB).