EMAT30008: Scientific Computing

Week 24: PDEs in 2D

Overview

In the final week of the unit, we will use finite differences to solve PDEs in two dimensions. Particular focus will be placed on solving Poisson's equation. For 2D problems, the number of unknowns and hence the size of the linear systems quickly grow as the number of grid points is increased. Hence, sparse matrices can lead to substantial reductions in the memory that needed to numerically solve the PDE.

Downloadable code

On Blackboard you can find a Python file that solves

2ux2+2uy2+q(x,y)=0\begin{aligned} \frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2} + q(x,y) = 0 \end{aligned}

on the rectangular domain axba \leq x \leq b, cydc \leq y \leq d with u=0u = 0 on all of the boundaries. This code uses NumPy to solve the linear system. The matrix D\mathbf{D} for the discretised Laplacian is stored as a dense NumPy array.

Exercises

The exercises below can be completed in any order.

  1. Write your own 2D Poisson solver that uses SciPy's root function.

  2. Extend the code on Blackboard so that non-homogeneous Dirichlet boundary conditions can be used. For example, replace the boundary condition at x=ax = a with u(a,y)=(yc)(yd)u(a, y) = (y-c)(y-d).

  3. Extend the code on Blackboard so that sparse matrices are used.

  4. Use memory profiling to determine how much memory is required to solve the 2D Poisson equation. How does the memory usage depend on the number of grid points?

  5. Write a solver for the 2D diffusion equation

    ut=D(2ux2+2uy2)\begin{aligned} \frac{\partial u}{\partial t} = D\left(\frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2}\right) \end{aligned}

    with u=0u = 0 at the boundaries and u(x,y,0)=1u(x,y,0) = 1. Beware, the time-step restrictions for 2D problems are different than for 1D problems!