09. Example Equations

In this section, we present a set of example equations to illustrate the concepts discussed in earlier sections. While the definitions.f90 file and settings.input configuration are provided in the solver/examples directory, it is strongly recommended that you attempt to recreate these results independently. This hands-on approach will enhance your understanding of the numerical methods and their implementation. We have included a few MATLAB codes in the directory solver if you want to visualise your solutions.

Example 1: Singular Perturbation Boundary Value Problem

The inputs for this example are found in this directory. Use the code to solve the equation

$$\epsilon \frac{\partial^2 u}{\partial x^2} - \frac{\partial u}{\partial x} = -1, \quad u(0) = 1,\quad u(1) = 3$$

See here for more details on this type of equation. This is a linear equation in one dimension. The general solution is given by $$u(x) = 1 + x + \frac{e^{x/\epsilon}-1}{e^{1/\epsilon}-1}.$$

• Without grid stretching, solve this equation in the cases of $\epsilon = 1,0.01$ and $0.001$ with $n_x=10$ as the number of grid points.
• What happens in the region past $x>0.8$ as $\epsilon$ gets smaller?
• Use grid stretching to cluster more points closer to end point $1$ of the domain in the cases of $\epsilon = 1,0.1,0.01$ and $0.001$. What can you conclude?

Example 2: Non-linear time dependent equation

The inputs for this example are found in this directory.

Use the code to solve the equation

$$\frac{\partial u}{\partial t} = \frac{\partial^2 u}{\partial x^2} + (1-u)u$$

on the domain $$0\leq x \leq 1$$

with boundary conditions

$$\frac{\partial u}{\partial x}(t,x=0) = 0$$ $$\frac{\partial u}{\partial x}(t,x=1) = 0$$

for the solution $u:=u(t,x)$. The initial conditions are given by

$$u(t=0,x) = \cos ^2(2 \pi x)$$

on time domain $[0, 1]$.

Things to note:

There are two interesting behaviours that occur from the solution: 1. First the solution becomes constant very rapidly. This is due to the diffusive term $\frac{\partial^2 u}{\partial x^2}$. 2. Secondly, the solution tends to $1$. This is because of the non-linear term $(1-u)u$: the time dependent equation is given by

$$\frac{\partial^2 u}{\partial x^2} + (1-u)u$$

and as $u$ becomes constant the term $$\frac{\partial^2 u}{\partial x^2}$$

becomes zero. We are left with

$$(1-u)u = 0$$

with the solution $u=1$ being an attractor i.e the solution $u=1$ is dominate over $u=0$.

Example 3: Predator-Prey Model

The inputs for this example are found in this directory.

Attempt to solve the following coupled equation

$$\epsilon (\frac{\partial^2 u}{\partial x^2}+\frac{\partial^2 u}{\partial y^2}) + u(1-u) - v\frac{u}{1+u} = \frac{\partial u}{\partial t}$$ $$(\frac{\partial^2 v}{\partial x^2}+\frac{\partial^2 v}{\partial y^2}) + v\frac{u}{1+u} = \frac{\partial v}{\partial t},$$

where $\epsilon = 0.05$, on the domain

$$0\leq x\leq 1$$

with boundary conditions

$$\frac{\partial u}{\partial x}(t,x=0) = 0, \quad \frac{\partial u}{\partial x}(t,x=1) = 0$$ $$\frac{\partial v}{\partial x}(t,x=0) = 0, \quad \frac{\partial v}{\partial x}(t,x=1) = 0$$

and time domain $$0\leq t \leq 0.5.$$

Use the initial conditions

$$u(t=0,x) = e^{-x^2}\cos^2(4\pi x)$$

$$v(t=0,x) = e^{(x+0.8)^2}(0.2 + 0.1\sin(10\pi x))^2$$

Exercises:

• How do the initial conditions change the behaviour of the predator and prey?
• Can you set up the system so that the prey die off?

Example 4: Turing Patterns

The inputs for this example are found in this directory.

Attempt to solve the following:

$$\epsilon (\frac{\partial^2 u}{\partial x^2}+\frac{\partial^2 u}{\partial y^2}) - u + u^2v = \frac{\partial u}{\partial t}\\(\frac{\partial^2 v}{\partial x^2}+\frac{\partial^2 v}{\partial y^2}) + \mu(1-u^2 v) = \frac{\partial v}{\partial t},$$

where $\epsilon = 0.05$ and $mu = 1.8$, on the domains

$$0\leq x\leq 20\\ \quad 0\leq y\leq 20$$

and time domain

$$0\leq t \leq 10.$$

Use the boundary conditions:

$$\frac{\partial u }{\partial x} (t,x=0,y) = \frac{\partial u }{\partial x}(t,x=20,y) = \frac{\partial u }{\partial y}(t,x,y=0) = \frac{\partial u }{\partial y}(t,x,y=20) = 0$$

$$\frac{\partial v }{\partial x}(t,x=0,y) = \frac{\partial v }{\partial x}(t,x=20,y) = \frac{\partial v }{\partial y}(t,x,y=0) = \frac{\partial v }{\partial y}(t,x,y=20) = 0$$

and initial conditions given by

$$u(t=0,x,y) = 1 + 0.01(2r-1)\\ \quad v(t=0,x,y) = 1$$

where $r$ is a random number between $0$ and $1$. Be careful with the number of spatial grid points; as we are in two dimensions the sizes of the matrices get huge very quickly and slow down the computation. This can be addressed with the use of parallel computing and openmp. We used $n_x=n_y=20$.

Hint: use the following intrinsic functions for r

  Call random_seed()
  Call random_number(r)

We note that for this particular problem, the small parameter $\epsilon$ that multiplies the second derivatives in the first equation does not cause a boundary layer to form in the PDE. Instead, the pattern made by $u$ will be less diffuse than the pattern made by $v$ all across the domain.

For a greater discussion see here.