```python
# Initialisation Cell
# You should always put imported modules here
from math import *
import numpy as np
import numpy.testing as nt
import numpy.linalg as LA
from matplotlib import pyplot as plt
import seaborn as sns
np.set_printoptions(suppress=True, precision=7)
```

# APPM4057 - November 2019

# Exam Office Info

| Examiner          | Name              |
|-------------------|-------------------|
| Internal Examiner | Prof B. A. Jacobs |
| External Examiner | Prof N. Hale (Stellenbosch) |

# CDEs Honours - Exam


## Instructions

* Read all the instructions carefully.
* The exam consists of **80 Marks**, with **three and a half** hours available.
* The written section is to be answered in the book provided.
* You must only access Moodle **TESTS** and NOT Moodle.
* The programming section is to be answered within this Jupyter notebook and resubmitted to Moodle **TESTS**.
* Do not rename the notebook, simply answer the questions and resubmit the file to Moodle.
* The moodle submission link will expire at exactly **12:30** and **NO** late submission will be accepted. Please make sure you submit timeously!
* The **Numpy** and **Matplotlib** documentation has been downloaded and is open on your current machine.
* **NB!!!** Anyone caught using Moodle (and its file storing), flash drives or external notes will be awarded zero and reported to the legal office.

## Written Section

* Answer the following questions in the answer book provided.

### Question 1 - 10 Marks

Discuss each of the below listed terms and their relation to one another (i.e. the theorem which describes this relationship):

* Consistency,
* Stability, 
* Convergence,

Be sure to give potential methods of analysing consistency, stability and convergence of a  linear PDE.

### Question 2 - 10 Marks

Analyse the stability of the difference scheme given by:
$$
u^{m}_{n} = -\beta u^{m+1}_{n-1} + (1+2\beta) u^{m+1}_{n} - \beta u^{m+1}_{n+1},
$$
using the discrete Fourier Transform:
$$
\hat{u}^{m+1}(\xi)=\frac{1}{\sqrt{2 \pi}} \sum_{n=-\infty}^{\infty} e^{-i n \xi} u_{n}^{m+1}.
$$

### Question 3 - 10 Marks

Consider the Crank-Nicolson (CN) scheme applied to the heat equation:
\begin{equation}
u_t = \nu u_{xx}.\label{eq:heat}
\end{equation}

(a) [3 Marks] Write out the CN scheme for the heat equation given by above.

(b) [7 Marks] Investigate the consistency of the scheme derived in Question 3(a).

### Question 4 - 10 Marks

Consider the following two-dimensional wave equation:
$$
\dfrac{\partial^2 u}{\partial t^2} = \nu^2 \left(\dfrac{\partial^2 u}{\partial x^2} + \dfrac{\partial^2 u}{\partial y^2}\right).\label{eq:wave}
$$

(a) [5 Marks] Set up the system of coupled first-order conservative partial differential equations by structuring new dependent variables $r$, $s$ and $l$, which can be used to solve the wave equation given above via the Lax-Wendroff finite difference scheme.

(b) [5 Marks] Provide the system in vector notation where the system is represented as:
$$
\dfrac{\partial \bf{U}}{\partial t} + \nabla \bf{F}(\bf{U}) = 0, 
$$
where $\bf{U} = [r,s,l]$.

## Programming Section 

### Question 1 - 15 Marks

(a) [10 Marks] Consider the wave equation given by:
$$
u_{tt} = c^2u_{xx},
$$
$$
u(x, 0) = x(1 - x), \ \ \ u_t(x, 0) = 0, \ \ \ u(0, t) = u(1, t) = 0.
$$
Write a function that implements a centered space, centered time scheme. Specifically, your function should take as inputs, the stepsizes `dx` and `dt`, the wave speed `c` (not squared valued), the number of iterations to perform `N` and the left and right endpoints `a`, `b`.

Your function should return the solution matrix, containing the wave profile at each time iteration, i.e. `N + 1` rows.


```python
def wave_eq(dx, dt, c, N, ic, a, b):
    # YOUR CODE HERE
    raise NotImplementedError()
    
```


```python
# Run this test cell to check your code
# Do not delete this cell
# 1 mark
# Unit test
dx = 0.1
c  = 0.2
tf = 8 
N  = 20
dt = (tf - 0)/N
a  = 0 
b  = 1
ic = lambda x: x*(1 - x)
tans = np.array([0.  , 0.09, 0.16, 0.21, 0.24, 0.25, 0.24, 0.21, 0.16, 0.09, 0.  ])
U = wave_eq(dx, dt, c, N, ic, a, b)
nt.assert_array_almost_equal(tans, U[0, :])
print('Test case passed!!!')
```


```python
# Run this test cell to check your code
# Do not delete this cell
# 1 mark
# Unit test
dx = 0.1
c  = 0.2
tf = 8 
N  = 20
dt = (tf - 0)/N
a  = 0 
b  = 1
ic = lambda x: x*(1 - x)
tans = np.array([0.    , 0.0836, 0.1536, 0.2036, 0.2336, 0.2436, 0.2336, 0.2036,
       0.1536, 0.0836, 0.    ])
U = wave_eq(dx, dt, c, N, ic, a, b)
nt.assert_array_almost_equal(tans, U[1, :])
print('Test case passed!!!')
```


```python
# Hidden test
# No output will be produced
# 8 marks
```

(b) [5 Marks] Continuing with Question 1(a), given the input below, produce the 3D surface plot, illustrating the evolution of the wave profile over time.


```python
dx    = 0.001
c     = 0.2
tf    = 8 
N     = 2000
dt    = (tf - 0)/N
a     = 0 
b     = 1
ic    = lambda x: x*(1 - x)
xvals = np.linspace(a, b, int((b - a)/dx) + 1)
U     = wave_eq(dx, dt, c, N, ic, a, b)
```


```python
# 5 Marks
# YOUR CODE HERE
raise NotImplementedError()
```

### Question 2 - 10 Marks


(a) [8 Marks] Write a function implements a finite difference scheme to solve the partial differential equation below. The function should implement an explicit scheme which is forward difference in time and central difference in space, specifically first order in time and second order in space. It should take as inputs a time step `dt`, a spatial step `dx`, a number of time-steps `N` to perform, the coefficients `D` and `nu`, an initial function `ic` passed as a handle and boundary values `alpha` and `beta`. The function should output the solution space matrix `u`, as well as the `xvals` (i.e. vector of $x$ steps) and `tvals` (i.e the vector of time steps). The PDE is given below:


$$v_t + \nu v_x = D v_{xx},$$

$$\text{BC}: v(0, t) = \alpha = v(1, t) = \beta = 0;\quad\text{IC}: v(x, 0) = \sin(8\pi x); \quad t > 0, \ x \in [0, 1]$$


```python
def heatEqn(dt, dx, N, ic, D, alpha, beta, nu):
    # YOUR CODE HERE
    raise NotImplementedError()
   

```


```python
# Run this test cell to check your code
# Do not delete this cell
# 1 mark
# Unit test
dx = 0.01
dt = 0.01
N  = 2
D  = 0.1
nu = 1
alpha = 0
beta  = 0
ic = lambda x: np.sin(8*np.pi*x) 
u, xx, tt = heatEqn(dt, dx, N, ic, D, alpha, beta, nu)
tans = np.array(  [ 0.       ,  0.2486899,  0.4817537,  0.6845471,  0.8443279,
                    0.9510565,  0.9980267,  0.9822873,  0.9048271,  0.7705132,
                    0.5877853,  0.3681246,  0.1253332, -0.1253332, -0.3681246,
                   -0.5877853, -0.7705132, -0.9048271, -0.9822873, -0.9980267,
                   -0.9510565, -0.8443279, -0.6845471, -0.4817537, -0.2486899,
                   -0.       ,  0.2486899,  0.4817537,  0.6845471,  0.8443279,
                    0.9510565,  0.9980267,  0.9822873,  0.9048271,  0.7705132,
                    0.5877853,  0.3681246,  0.1253332, -0.1253332, -0.3681246,
                   -0.5877853, -0.7705132, -0.9048271, -0.9822873, -0.9980267,
                   -0.9510565, -0.8443279, -0.6845471, -0.4817537, -0.2486899,
                   -0.       ,  0.2486899,  0.4817537,  0.6845471,  0.8443279,
                    0.9510565,  0.9980267,  0.9822873,  0.9048271,  0.7705132,
                    0.5877853,  0.3681246,  0.1253332, -0.1253332, -0.3681246,
                   -0.5877853, -0.7705132, -0.9048271, -0.9822873, -0.9980267,
                   -0.9510565, -0.8443279, -0.6845471, -0.4817537, -0.2486899,
                   -0.       ,  0.2486899,  0.4817537,  0.6845471,  0.8443279,
                    0.9510565,  0.9980267,  0.9822873,  0.9048271,  0.7705132,
                    0.5877853,  0.3681246,  0.1253332, -0.1253332, -0.3681246,
                   -0.5877853, -0.7705132, -0.9048271, -0.9822873, -0.9980267,
                   -0.9510565, -0.8443279, -0.6845471, -0.4817537, -0.2486899,
                    0.       ])
nt.assert_array_almost_equal(tans, u[0])
print('Test case passed!!!')
```


```python
# Run this test cell to check your code
# Do not delete this cell
# 1 mark
# Unit test
dx = 0.01
dt = 0.01
N  = 2
D  = 0.1
nu = 1
alpha = 0
beta  = 0
ic = lambda x: np.sin(8*np.pi*x) 
u, xx, tt = heatEqn(dt, dx, N, ic, D, alpha, beta, nu)
tans = np.array(  [ 0.       ,  2.4511655, -0.1252406, -0.0825335, -0.0346406,
                    0.015429 ,  0.064529 ,  0.1095745,  0.1477351,  0.1766128,
                    0.1943934,  0.1999595,  0.1929614,  0.1738388,  0.1437933,
                    0.1047128,  0.0590527,  0.0096822, -0.0402967, -0.0877437,
                   -0.1296773, -0.1634629, -0.1869775, -0.1987436, -0.1980219,
                   -0.1848578, -0.1600783, -0.1252406, -0.0825335, -0.0346406,
                    0.015429 ,  0.064529 ,  0.1095745,  0.1477351,  0.1766128,
                    0.1943934,  0.1999595,  0.1929614,  0.1738388,  0.1437933,
                    0.1047128,  0.0590527,  0.0096822, -0.0402967, -0.0877437,
                   -0.1296773, -0.1634629, -0.1869775, -0.1987436, -0.1980219,
                   -0.1848578, -0.1600783, -0.1252406, -0.0825335, -0.0346406,
                    0.015429 ,  0.064529 ,  0.1095745,  0.1477351,  0.1766128,
                    0.1943934,  0.1999595,  0.1929614,  0.1738388,  0.1437933,
                    0.1047128,  0.0590527,  0.0096822, -0.0402967, -0.0877437,
                   -0.1296773, -0.1634629, -0.1869775, -0.1987436, -0.1980219,
                   -0.1848578, -0.1600783, -0.1252406, -0.0825335, -0.0346406,
                    0.015429 ,  0.064529 ,  0.1095745,  0.1477351,  0.1766128,
                    0.1943934,  0.1999595,  0.1929614,  0.1738388,  0.1437933,
                    0.1047128,  0.0590527,  0.0096822, -0.0402967, -0.0877437,
                   -0.1296773, -0.1634629, -0.1869775, -0.1987436,  2.164532 ,
                    0.       ])
nt.assert_array_almost_equal(tans, u[-1])
print('Test case passed!!!')
```


```python
# Hidden test
# No output will be produced
# 6 marks
```

(b) [2 Marks] Continuing on from Question 3(a). Plot the profile of each time step on the same set of axes, given the following inputs (be sure to use a legend indicating which iteration and curve is which):


```python
dx    = 0.01
dt    = 0.01
N     = 2
D     = 0.1
nu    = 1
alpha = 0
beta  = 0
ic    = lambda x: np.sin(8*np.pi*x) 
```


```python
# 2 Marks
# YOUR CODE HERE
raise NotImplementedError()
```

### Question 3 - 15 Marks

(a) [10 Marks] Consider the PDE given by:
$$
v_t = D\left( v_{xx} + v_{yy}\right), 
$$
$$
v(-1, y, t) = v(1, y, t) = 0,\ v(x, y, 0) = 0,  
$$
\begin{equation}
v(x, -1, 0) = v(x, 1, 0) = 100, \ \ \ x \in (-1, 1)\ \ \text{i.e. endpoints excluded}
\end{equation}

Periodic boundary conditions should be imposed on $y = -1, y = 1$, i.e.$v(x,-1,t)=v(x,1,t)$ for all $t$ and $x$.


Write a function which implements an ADI scheme to solve the PDE given above. Specifically, the function should take as inputs, the step-sizes, `dx` and `dy`, the time-step `dt`, the constant `D`, and the number of iterations to perform `N`.

The function should return the 2D heatmap for each time-step, that is, a 3D array of `N + 1` time-steps.


```python
def heat2D_ADI(dx, dy, dt, D, N):
    # YOUR CODE HERE
    raise NotImplementedError()
```


```python
# Run this test cell to check your code
# Do not delete this cell
# 1 mark
# Unit test
dx   = 0.2
dy   = 0.2
dt   = 0.1
D    = 0.1
N    = 5
tans = np.array([[  0., 100., 100., 100., 100., 100., 100., 100., 100., 100.,   0.],
       [  0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.],
       [  0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.],
       [  0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.],
       [  0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.],
       [  0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.],
       [  0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.],
       [  0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.],
       [  0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.],
       [  0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.],
       [  0., 100., 100., 100., 100., 100., 100., 100., 100., 100.,   0.]])
U    = heat2D_ADI(dx, dy, dt, D, N)
nt.assert_array_almost_equal(tans, U[:, :, 0])
print('Test case passed!!!')
```


```python
# Run this test cell to check your code
# Do not delete this cell
# 1 mark
# Unit test
dx   = 0.2
dy   = 0.2
dt   = 0.1
D    = 0.1
N    = 5
tans = np.array([[ 0.       , 65.3025091, 78.3317661, 79.6479746, 79.780803 ,
        79.7928783, 79.780803 , 79.6479746, 78.3317661, 65.3025091,
         0.       ],
       [ 0.       , 14.8641016, 17.82981  , 18.1294042, 18.1596385,
        18.1623871, 18.1596385, 18.1294042, 17.82981  , 14.8641016,
         0.       ],
       [ 0.       ,  1.5015808,  1.8011785,  1.8314437,  1.834498 ,
         1.8347756,  1.834498 ,  1.8314437,  1.8011785,  1.5015808,
         0.       ],
       [ 0.       ,  0.1517061,  0.1819747,  0.1850325,  0.185341 ,
         0.1853691,  0.185341 ,  0.1850325,  0.1819747,  0.1517061,
         0.       ],
       [ 0.       ,  0.0154802,  0.0185689,  0.0188809,  0.0189124,
         0.0189152,  0.0189124,  0.0188809,  0.0185689,  0.0154802,
         0.       ],
       [ 0.       ,  0.003096 ,  0.0037138,  0.0037762,  0.0037825,
         0.003783 ,  0.0037825,  0.0037762,  0.0037138,  0.003096 ,
         0.       ],
       [ 0.       ,  0.0154802,  0.0185689,  0.0188809,  0.0189124,
         0.0189152,  0.0189124,  0.0188809,  0.0185689,  0.0154802,
         0.       ],
       [ 0.       ,  0.1517061,  0.1819747,  0.1850325,  0.185341 ,
         0.1853691,  0.185341 ,  0.1850325,  0.1819747,  0.1517061,
         0.       ],
       [ 0.       ,  1.5015808,  1.8011785,  1.8314437,  1.834498 ,
         1.8347756,  1.834498 ,  1.8314437,  1.8011785,  1.5015808,
         0.       ],
       [ 0.       , 14.8641016, 17.82981  , 18.1294042, 18.1596385,
        18.1623871, 18.1596385, 18.1294042, 17.82981  , 14.8641016,
         0.       ],
       [ 0.       , 65.3025091, 78.3317661, 79.6479746, 79.780803 ,
        79.7928783, 79.780803 , 79.6479746, 78.3317661, 65.3025091,
         0.       ]])

U = heat2D_ADI(dx, dy, dt, D, N)
nt.assert_array_almost_equal(tans, U[:, :, 1])
print('Test case passed!!!')
```


```python
# Hidden test
# No output will be produced
# 8 marks
```

(b) [2 Marks] Plot a heatmap of the final time-step using the ADI function in Question 2(a), given the input below:


```python
dx = 0.1
dy = 0.1
dt = 0.1
D  = 0.2
N  = 5
U  = heat2D_ADI(dx, dy, dt, D, N) 
```


```python
# 2 Marks
# YOUR CODE HERE
raise NotImplementedError()
```

(c) [3 Marks] Modify your function from Question 3(a) to plot the two tridiagonal ADI matrices $A$ and $B$.


```python
#3 Marks
# YOUR CODE HERE
raise NotImplementedError()
```

a518477516ad7f55a08d29fe10221da7275f8b63
ipynb
Jupyter Notebook
2019/release/Exam/exam.ipynb
mattw89/APPM4057-Computational-Differential-Equations
dd588de401c9f1b26f13ed1cf7a85fec5a891f9f
CC-BY-4.0
2019/release/Exam/exam.ipynb
mattw89/APPM4057-Computational-Differential-Equations
dd588de401c9f1b26f13ed1cf7a85fec5a891f9f
CC-BY-4.0
2019/release/Exam/exam.ipynb
mattw89/APPM4057-Computational-Differential-Equations
dd588de401c9f1b26f13ed1cf7a85fec5a891f9f
CC-BY-4.0
Qwen/Qwen-72B
1. YES
2. YES