Assignment 3

Authors: Joaquín Pelle, Oscar Reula, Pablo Montes

Click here to download this notebook

Important 1:

This notebook is designed to run on Julia 1.6 or later.

Important 2:

Rename the file as: name_assignment_3.ipynb

All generated code and the presentation must be included in this notebook.

The objectives of this assignment are as follows:

  1. To understand wave propagation and relate it to the need to impose boundary conditions according to the direction of propagation.

  2. To comprehend the scheme of imposing boundary conditions using finite difference operators with the Summation By Parts property.


1) In the directory, you will find a Julia notebook that allows you to evolve the advection equation by providing a boundary condition at the edge where the solution enters the integration region. Try different initial and boundary conditions. What happens if you change the sign of the velocity?

2) Note that the convergence of the method, in the $L^2$ norm, is of second order, even though the approximation at the boundary is first order.

3) The fourth order Summation By Parts finite difference operator has the following matrix entries near $1$ and $N$:

using LinearAlgebra
Head = [
   -1.4117647059   1.7352941176  -0.2352941176  -0.0882352941   0.0           0.0
   -0.5            0.0            0.5            0.0            0.0           0.0
    0.0930232558  -0.6860465116   0.0            0.6860465116  -0.0930232558  0.0
    0.0306122449   0.0           -0.6020408163   0.0            0.6530612245 -0.0816326531
]
4×6 Matrix{Float64}:
 -1.41176     1.73529   -0.235294  -0.0882353   0.0         0.0
 -0.5         0.0        0.5        0.0         0.0         0.0
  0.0930233  -0.686047   0.0        0.686047   -0.0930233   0.0
  0.0306122   0.0       -0.602041   0.0         0.653061   -0.0816327

Plus a standard band matrix part with bandwidth 5. Note that the bottom part of the matrix is not the same as the Head, but instead corresponds to an appropriate transposition and sign change of it.

middle = [0.0833333333333 -0.6666666666666 0  0.6666666666666 -0.0833333333333]
1×5 Matrix{Float64}:
 0.0833333  -0.666667  0.0  0.666667  -0.0833333

On the other hand, the inner product is given by:

function h4(N)
    h = ones(N)
    # Set specialized values at the boundaries
    h[1]   = 17/48
    h[2]   = 59/48
    h[3]   = 43/48
    h[4]   = 49/48
    h[end]     = h[1]
    h[end-1]   = h[2]
    h[end-2]   = h[3]
    h[end-3]   = h[4]
    return Diagonal(h)
end

h4(10)
10×10 Diagonal{Float64, Vector{Float64}}:
 0.354167   ⋅        ⋅         ⋅        ⋅   …   ⋅         ⋅        ⋅ 
  ⋅        1.22917   ⋅         ⋅        ⋅       ⋅         ⋅        ⋅ 
  ⋅         ⋅       0.895833   ⋅        ⋅       ⋅         ⋅        ⋅ 
  ⋅         ⋅        ⋅        1.02083   ⋅       ⋅         ⋅        ⋅ 
  ⋅         ⋅        ⋅         ⋅       1.0      ⋅         ⋅        ⋅ 
  ⋅         ⋅        ⋅         ⋅        ⋅   …   ⋅         ⋅        ⋅ 
  ⋅         ⋅        ⋅         ⋅        ⋅       ⋅         ⋅        ⋅ 
  ⋅         ⋅        ⋅         ⋅        ⋅      0.895833   ⋅        ⋅ 
  ⋅         ⋅        ⋅         ⋅        ⋅       ⋅        1.22917   ⋅ 
  ⋅         ⋅        ⋅         ⋅        ⋅       ⋅         ⋅       0.354167

Construct the fourth-order operator and verify that it converges with order 4 in the $L^2$ norm (even though near the boundary the convergence drops to second order). Also, verify that it satisfies the Summation By Parts (SBP) property.

4) Now consider the wave equation, $\phi_{tt} = c^2 \phi_{xx}$. Decompose it into a first-order system with one mode propagating to the left and one to the right. Hint: Determine what equation is satisfied by the combination $V_+ := \phi_t - c\phi_x$. Use this decomposition to specify boundary conditions so that the problem is well-posed (but not overdetermined). Consider the following cases:

  1. Homogeneous boundary conditions: nothing enters, but the initial data is nonzero.

  2. Reflecting boundary conditions: whatever arrives as one mode enters as the other. To impose this condition, choose the combination of modes at the boundary so that $\phi$ remains constant in time.

5) For the adventurous: Consider two adjacent regions, $L_L$ and $L_R$, in which the wave equation evolves but with different speeds $(c_L, c_R)$ (imagine that the regions are strings made of different materials or densities). Treat the interface between the media as a place where boundary conditions are specified according to what comes from the other region. Choose these conditions so that $\phi_t$ and $\phi_x$ are continuous across the interface. Implement code that reflects the situation described.

  • Implement boundary conditions on the left side of $L_R$ so that a sinusoidal wave enters.

  • Implement boundary conditions on the right side of $L_L$ so that the outgoing wave is reflected and becomes incoming.

Click here to see an example gif of what you should achieve (the scale of the waves is off, focus on there being a reflection and a transmission):