Solving the Poisson equation using Python

This  post is part of the CFDPython series that shows how to solve the Navier Stokes equations with finite difference method by use of Python. This is the last step to the small solver we want to create. The Poisson equation is actually the Laplace equation to which we add a source term to the right hand side:

∂2p∂x2+∂2p∂y2=b

2D Poisson equation solve for 4 source terms

To better understand how this appears, we must return at the Navier Stokes equations; for a incompressible fluid we must resolve the continuity and the momentum equations; the continuity equation creates a kinematic constraint which implies that div(U) should be zero everywhere in the domain and is not coupled with the pressure; that is why we cannot use the continuity directly in the developing a solving procedure.
Actually the continuity is implicitly assured by creating a pressure field from the momentum equations; the equations we will solve here is a simplification of the actual Poisson equation; the full equation can be seen here.
Coming back to our problem we see we can use the step 9 solution with some alterations; the source term will give a value on the field inside the domain; the diffusion will try to smooth this value.
The discretization of the equation is done using the central difference scheme:

pni+1,j−2pni,j+pni−1,jΔx2+pni,j+1−2pni,j+pni,j−1Δy2=bni,j

And the computing scheme is the following:

pni,j=(pni+1,j+pni−1,j)Δy2+(pni,j+1+pni,j−1)Δx2−bni,jΔx2Δy22(Δx2+Δy2)

For this equation we must define initial conditions, boundary condition and source terms:
- IC: pressure is zero everywhere
- BC: pressure is zero at the domain boundaries ( p = 0 at x = 0, 2 and y = 0, 2)
- the source term consist in two initial spikes inside the domain:
bi,j=100 at i = nx/4 and j = ny/4
bi,j=−100 at i = 3*nx/4 and j = 3*ny/4
bi,j=0 everywhere else.
The Python code for this problem is given below:

And the solution after 1000 steps (pseudo-time) is the following:

At this moment all the components from the Navier Stokes equations have been solved by the use of Python. All codes may be united to create a 2D finite difference solver. This will actually be used next to solve some basic problems of fluid dynamics: the lid driven cavity flow and the viscous flow in a pipe.
The codes shown here are available on my GitHub page; if there are any questions please leave a comment below.