{
"cells": [
{
"cell_type": "markdown",
"id": "81b1d2f0",
"metadata": {},
"source": [
"---\n",
"title: \"Wave Equation\"\n",
"subtitle: \"Lecture 20\"\n",
"date: 2026-04-22\n",
"abstract-title: Overview\n",
"format:\n",
" html:\n",
" other-links:\n",
" - text: This notebook\n",
" href: L20.ipynb\n",
"---"
]
},
{
"cell_type": "markdown",
"id": "63cf0867",
"metadata": {},
"source": [
"::: {.callout-note}\n",
"\n",
"I encourage you to play around with the juptyer notebook for this lecture - you can copy the code with the ```this notebook``` button on the side of this page.\n",
"\n",
":::\n",
"\n",
"::: {.callout-warning}\n",
"\n",
"These notes are mainly a record of what we discussed and are not a substitute for attending the lectures and reading books! If anything is unclear/wrong, let me know and I will update the notes.\n",
"\n",
":::\n",
"\n",
"::: {.callout-tip}\n",
"\n",
"This chapter is mostly based on \n",
"\n",
"* @Burden Chapter 11, \n",
"* @FNC Chapter 10, \n",
"* Lecture 16: @FNC [Galerkin method](https://fncbook.com/galerkin/), \n",
"* Lecture 17: Chapter 14.4 of @Suli: Theorem 14.6. \n",
"* Lecture 19: @Olver2014 Chapter 5\n",
"\n",
":::"
]
},
{
"cell_type": "code",
"execution_count": 25,
"id": "f8408d86",
"metadata": {},
"outputs": [],
"source": [
"#| echo: false\n",
"\n",
"using Plots, LaTeXStrings, LinearAlgebra, Interpolations, SparseArrays\n",
"\n",
"# IF YOU ARE READING THIS, THANK YOU\n",
"# FOR DOWNLOADING THE LECTURE NOTES.\n",
"# THE FIRST PERSON TO LET ME KNOW,\n",
"# WINS SOME CHOCOLATES OF YOUR CHOICE"
]
},
{
"cell_type": "markdown",
"id": "81473756",
"metadata": {},
"source": [
"## Recap\n",
"\n",
"*Assignment 6:* You are looking at simple finite difference approximations to the linear transport equation $u_t + c u_x = 0$. You will see convergence under an assumption on the *Courant number* $\\sigma := c\\frac{\\tau}{h}$.\n",
"\n",
"*Last time:* Convergence of FTCS and BTCS finite difference methods for solving the heat equation $u_t - \\alpha^2 u_{xx} = 0$. For FTCS, we saw that we have convergence when $0 \\leq \\mu \\leq \\frac{1}{2}$ where $\\mu := \\tau (\\frac{\\alpha}{h})^2$ is the *diffusion number* (and we get faster convergence when $\\mu = \\frac16$). Then, we saw that BTCS is unconditionally stable.\n",
"\n",
"Today, we look at the wave equation: \n",
"\n",
"\\begin{align}\n",
" &\\square u := u_{tt} - c^2 u_{xx} = 0 \\tag{wave}\\\\\n",
" &\\begin{split}\n",
" &u(x,0) = u_0(x) \\\\\n",
" &u_t(x,0) = v_0(x)\n",
" \\end{split} \\tag{IC} \\\\\n",
" &\\begin{split}\n",
" &u(a,t) = g_a(t) \\\\\n",
" &u(b,t) = g_b(t)\n",
" \\end{split} \\tag{BC}\n",
"\\end{align}\n",
"\n",
"where $c$ is the wave speed. (Note: because the equation is second order in time, we require two initial conditions).\n",
"\n",
"::: {#exr-}\n",
"\n",
"(Ignoring the initial and boundary conditions) show that $u_0(x \\pm c t)$ solve the $\\text{(wave)}$. That is, the wave equation supports waves travelling in both directions.\n",
"\n",
":::"
]
},
{
"cell_type": "markdown",
"id": "793aabd3",
"metadata": {},
"source": [
"## Finite differences \n",
"\n",
"We again consider a uniform mesh:\n",
"\n",
"\\begin{align}\n",
" x_i &= a + \\frac{b-a}{m} i =: a + ih &&\\text{for }i=0,\\dots,m\\nonumber\\\\\n",
" t_j &= \\frac{T}{n} j =: j\\tau &&\\text{for }j=0,\\dots,n\n",
"\\end{align}\n",
"\n",
"and aim to approximate $u(x_i,t_j) \\approx u_{ij}$.\n",
"\n",
"The simplest approximation we can think of is to replace the second derivatives with our usual centred difference:\n",
"\n",
"\\begin{align}\n",
" \\frac\n",
" {u_{i,j-1} - 2 u_{ij} + u_{i,j+1}}\n",
" {\\tau^2}\n",
" %\n",
" - c^2 \\frac\n",
" {u_{i-1,j} - 2 u_{ij} + u_{i+1,j}}\n",
" {h^2}\n",
" %\n",
" = 0\n",
"\\end{align}\n",
"\n",
"which is equivalent to \n",
"\n",
"\\begin{align}\n",
" u_{i,j+1}\n",
" %\n",
" = \\sigma^2 u_{i-1,j} + 2( 1 - \\sigma^2) u_{ij} + \\sigma^2 u_{i+1,j} \n",
" %\n",
" - u_{i,j-1}\n",
"\\end{align}\n",
"\n",
"where $\\sigma := c\\frac{\\tau}{h}$. \n",
"\n",
"We can impose the boundary conditions by setting $u_{0,j} := g_a(t_j)$ and $u_{m,j} := g_b(t_j)$.\n",
"\n",
"::: {#exr-}\n",
"\n",
"Write the above as a matrix equation for the numerical solution in the interior mesh points $\\{ u_{i,j+1} \\}_{i=1}^{m-1}$.\n",
"\n",
":::\n",
"\n",
"There is a slight problem here! In order to compute $u_{i,j+1}$ we require the solution at the previous *two* time steps. In order to start the iteration we therefore need to have $\\{u_{i,0}\\}$ and $\\{u_{i,1}\\}$. As we have done so far, we can set $u_{i,0} := u_0(x_i)$ but we also need to specify $\\{u_{i,1}\\}$. We do this by using the other initial condition: $u_t(x,0) = v_0(x)$.\n",
"\n",
"Let $u$ be the exact solution to the wave equation. Then, we have \n",
"\n",
"\\begin{align}\n",
" u(x_i,t_1 ) &= u(x_i,\\tau) \n",
" %\n",
" = u(x_i,0) + \\tau u_t(x_i,0) + \\frac{\\tau^2}{2} u_{tt}(x_i,0) + O(\\tau^3) \\nonumber\\\\\n",
" %\n",
" &= u_0(x_i) + \\tau v_0(x_i) + \\frac{c^2 \\tau^2}{2} u_0''(x_i) + O(\\tau^3) \\nonumber\\\\\n",
" %\n",
" &= u_0(x_i) + \\tau v_0(x_i) + \\tfrac{\\sigma^2}{2} \\left[ u_0(x_{i-1}) - 2u_{0}(x_i) + u_0(x_{i+1}) \\right] \\nonumber\\\\\n",
" &\\qquad + O(\\tau^3) + O(h^2)\n",
"\\end{align} \n",
"\n",
"In particular, we could define $\\{u_{i,1}\\}$ by \n",
"\n",
"\\begin{align}\n",
" u_{i,1} := \\tfrac{\\sigma^2}{2} u_{i-1,0} +(1- \\sigma^2) u_{i,0} + \\tfrac{\\sigma^2}{2} u_{i+1,0} + \\tau v_0(x_i)\n",
"\\end{align}\n",
"\n",
"We implement this here:"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "905ef000",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"σ = 0.9\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"\u001b[36m\u001b[1m[ \u001b[22m\u001b[39m\u001b[36m\u001b[1mInfo: \u001b[22m\u001b[39mSaved animation to c:\\Users\\math5\\Math 5485\\Math5486\\pics\\wave-1.gif\n"
]
},
{
"data": {
"text/html": [
"
"
],
"text/plain": [
"Plots.AnimatedGif(\"c:\\\\Users\\\\math5\\\\Math 5485\\\\Math5486\\\\pics\\\\wave-1.gif\")"
]
},
"execution_count": 26,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# | code-fold: show\n",
"\n",
"function fdm_wave( m, n, c, u0, v0; a=0., b=1., f=0., g_a = 0., g_b = 0., T=.2 )\n",
" # space discretisation\n",
" h = (b-a)/m\n",
" x = collect(a:h:b)\n",
"\n",
" # time discretisation\n",
" τ = T/n \n",
" t = collect(0:τ:T)\n",
"\n",
" σ = c*τ/h\n",
" println( \"σ = $σ\" )\n",
" \n",
" u = zeros( m+1, n+1 )\n",
"\n",
" # BC\n",
" u[1,:] = typeof(g_a) <: Function ? g_a.(t) : fill(g_a,n+1)\n",
" u[m+1,:] = typeof(g_b) <: Function ? g_b.(t) : fill(g_b,n+1)\n",
"\n",
" # IC\n",
" u[:,1] = u0.(x)\n",
" u[2:m,2] = (σ^2/2)u0.(x[1:m-1]) + (1-σ^2)u0.(x[2:m]) +(σ^2/2)u0.(x[3:m+1]) + τ*v0.(x[2:m])\n",
"\n",
" for j ∈ 2:n\n",
" for i ∈ 2:m\n",
" source = typeof(f) <: Function ? f(x[i],t[i]) : f\n",
" u[i,j+1] = (σ^2)u[i-1,j] + 2(1-σ^2)u[i,j] + (σ^2)u[i+1,j] - u[i,j-1] + τ^2 * source\n",
" end\n",
" end\n",
"\n",
" return x,t,u\n",
"end\n",
"\n",
"c = 1.\n",
"\n",
"# initial condition\n",
"u0(x) = exp( -400(x-.3)^2 )\n",
"v0(x) = 0.\n",
"\n",
"# exact solution\n",
"g(z) = mod(z+1,2)-1 < 0 ? -u0(-(mod(z+1,2)-1)) : u0(mod(z+1,2)-1)\n",
"u_e(x,t) = (1/2)*(g(x-t)+g(x+t))\n",
"\n",
"x,t,u = fdm_wave( 90, 1000, c, u0, v0; T=10. )\n",
"\n",
"anim = @animate for j ∈ 1:5:length(t)\n",
" plot(x, u[:, j]; \n",
" title = L\"Wave equation $u_{tt} - u_{xx} = 0$\",\n",
" label = \"t = $(t[j])\", legend=:topright,\n",
" ylims = [-1,1], m=2)\n",
" plot!( x->u_e(x,t[j]); \n",
" label=\"exact solution\", linestyle=:dot)\n",
"end\n",
"gif(anim, \"pics/wave-1.gif\", fps=5)"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "42eb9cda",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"σ = 1.01\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"\u001b[36m\u001b[1m[ \u001b[22m\u001b[39m\u001b[36m\u001b[1mInfo: \u001b[22m\u001b[39mSaved animation to c:\\Users\\math5\\Math 5485\\Math5486\\pics\\wave-2.gif\n"
]
},
{
"data": {
"text/html": [
"
"
],
"text/plain": [
"Plots.AnimatedGif(\"c:\\\\Users\\\\math5\\\\Math 5485\\\\Math5486\\\\pics\\\\wave-2.gif\")"
]
},
"execution_count": 27,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# | code-fold: show\n",
"\n",
"x,t,u = fdm_wave( 101, 250, c, u0, v0; T=2.5 )\n",
"\n",
"anim = @animate for j ∈ 1:5:length(t)\n",
" plot(x, u[:, j]; \n",
" title = L\"Wave equation $u_{tt} - u_{xx} = 0$\",\n",
" label = \"t = $(t[j])\", legend=:topright,\n",
" ylims = [-1,1], m=2)\n",
" plot!( x->u_e(x,t[j]); \n",
" label=\"exact solution\", linestyle=:dot)\n",
"end\n",
"gif(anim, \"pics/wave-2.gif\", fps=10)"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "a3cec3a4",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"σ = 0.9\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"\u001b[36m\u001b[1m[ \u001b[22m\u001b[39m\u001b[36m\u001b[1mInfo: \u001b[22m\u001b[39mSaved animation to c:\\Users\\math5\\Math 5485\\Math5486\\pics\\wave-3.gif\n"
]
},
{
"data": {
"text/html": [
"
"
],
"text/plain": [
"Plots.AnimatedGif(\"c:\\\\Users\\\\math5\\\\Math 5485\\\\Math5486\\\\pics\\\\wave-3.gif\")"
]
},
"execution_count": 28,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# | code-fold: show\n",
"\n",
"c = 1.\n",
"\n",
"# initial condition\n",
"u0(x) = x*sin(2π*x) #abs(x-.5)<.75 ? (abs(x-.5)<.25 ? 0 : 1/2 ) : 0\n",
"v0(x) = 0.\n",
"\n",
"# exact solution\n",
"g(z) = mod(z+1,2)-1 < 0 ? -u0(-(mod(z+1,2)-1)) : u0(mod(z+1,2)-1)\n",
"u_e(x,t) = (1/2)*(g(x-t)+g(x+t))\n",
"\n",
"x,t,u = fdm_wave( 90, 1000, c, u0, v0; T=10. )\n",
"\n",
"anim = @animate for j ∈ 1:5:length(t)\n",
" plot(x, u[:, j]; \n",
" title = L\"Wave equation $u_{tt} - u_{xx} = 0$\",\n",
" label = \"t = $(t[j])\", legend=:topright,\n",
" ylims = [-1,1], m=2)\n",
" plot!( x->u_e(x,t[j]); \n",
" label=\"exact solution\", linestyle=:dot)\n",
"end\n",
"gif(anim, \"pics/wave-3.gif\", fps=5)"
]
},
{
"cell_type": "markdown",
"id": "e4da2d22",
"metadata": {},
"source": [
"# Courant–Friedrichs–Lewy (CFL) condition\n",
"\n",
"::: {#exr-}\n",
"\n",
"Use the change of variables $\\xi(x,t) := x-ct$ and $\\eta(x,t) := x+ct$ to rewrite the wave equation $\\square u = 0$ as $u_{\\xi\\eta} = 0$. Hence, prove the d'Alembert formula:\n",
"\n",
"\\begin{align}\n",
" u(x,t) = \\frac{u_0(x-ct) + u_0(x+ct)}{2}\n",
" %\n",
" +\\frac{1}{2c}\\int_{x-ct}^{x+ct} v_0(y)\\mathrm{d}y.\n",
"\\end{align}\n",
"\n",
":::\n",
"\n",
"The (exact) domain of dependence at $(x,t)$ is the set of all $y$ for which changing the initial data $u_0(y)$, $v_0(y)$ can affect the solution at $(x,t)$. \n",
"\n",
"::: {#exr-}\n",
"\n",
"Use d'Alembert formula to write down the domain of dpendence for the wave equation. \n",
"\n",
":::\n",
"\n",
"The numerical domain of dependence at $(x_i,t_j)$ is the set of all $x_k$ for which $u_{ij}$ is affected by $u_{k0}$.\n",
"\n",
"::: {#exr-}\n",
"\n",
"Consider the finite difference scheme as above: $u_{i,0} := u_0(x_i)$ and \n",
"\n",
"\\begin{align}\n",
" u_{i,1} \n",
" %\n",
" &= \\tau v_0(x_i) + \\tfrac{\\sigma^2}{2} u_{i-1,0} +(1- \\sigma^2) u_{i,0} + \\tfrac{\\sigma^2}{2} u_{i+1,0} \\nonumber\\\\\n",
" %\n",
" u_{i,j+1}\n",
" %\n",
" &= \\sigma^2 u_{i-1,j} + 2( 1 - \\sigma^2) u_{ij} + \\sigma^2 u_{i+1,j} \n",
" %\n",
" - u_{i,j-1}\n",
"\\end{align}\n",
"\n",
"where $\\sigma := c\\frac{\\tau}{h}$. What is the numerical domain of dependence at $(x_i,t_j)$ of this method?\n",
"\n",
":::\n",
"\n",
"The Courant–Friedrichs–Lewy (CFL) condition states that, for a numerical method to converge, the exact domain of dependence must be contained in the numerical domain of dependence as $h, \\tau \\to 0$.\n",
"\n",
"::: {#exr-}\n",
"\n",
"Under what condition on $\\sigma$ does the numerical scheme satisfy the CFL condition?\n",
"\n",
":::\n",
"\n",
"::: {#exr-}\n",
"\n",
"Does the CLF condition explain the instability of the method that we noticed above?\n",
"\n",
":::"
]
},
{
"cell_type": "markdown",
"id": "8f8141aa",
"metadata": {},
"source": [
"::: {#exr-}\n",
"\n",
"Complete the von Neumann stability analysis to determine when the numerical scheme is stable. \n",
"\n",
"Hint: Take $u_{i,j-1} = e^{\\mathrm{i}kx_i}$, $u_{i,j} =\\lambda e^{\\mathrm{i}kx_i}$, and $u_{i,j+1} =\\lambda^2 e^{\\mathrm{i}kx_i}$ and show that $\\lambda$ satisfies a quadratic equation giving $\\lambda = \\alpha \\pm \\sqrt{\\alpha^2 - 1}$. Show that when $\\sigma$ satisfies the CFL condition then $|\\lambda|=1$, and, when $\\sigma$ violates the CFL condition then one of the $\\lambda's$ has absolute value $> 1$.\n",
"\n",
"In this case the CFL condition distinguishes between stable and unstable constants $\\sigma$. This is not true in general - see A6.\n",
"\n",
":::"
]
},
{
"cell_type": "markdown",
"id": "81dc004c",
"metadata": {},
"source": [
"## Semi-discretisation\n",
"\n",
"In this lecture we have treated time in the same way as space. Another approach is the so called *method of lines*:\n",
"\n",
"Another general approach is semi-discretisation: we again take $x_i := a + \\frac{b-a}{m} i$ and $h := \\frac{b-a}{m}$ and now approximate $\\{u(x_i, t)\\}_{i=0}^m$ with the (time-dependent) vector $\\bm u(t)$. (This is known as the semi-discretisation step because we have discretised in space but not time). \n",
"\n",
"Then, if one is interested in solving the heat equation ($u_t = \\alpha^2 u_{xx}$), we can replace this with \n",
"\n",
"\\begin{align}\n",
" \\bm u'(t) &= \\alpha^2 D\\bm u(t)\n",
"\\end{align}\n",
"\n",
"where $D$ is a differentiation matrix approximating the second derivative in space. This is a linear, constant coefficient system of ODEs that we can solve directly using techniques from Chapter 1 (also see A2). \n",
"\n",
"::: {#exr-}\n",
"\n",
"Approximate the heat equation $\\square u := u_{tt} - c^2 u_{xx} = 0$ with a linear system of ODEs using the method of lines. \n",
"\n",
"Hint: you could define *(i)* $u_t = y$ and $y_t = c^2 u_{xx}$ or *(ii)* $u_t = z_x$ and $z_t = c^2 u_{x}$. (Here, *(ii)* gives *Maxwell's equations*).\n",
"\n",
":::"
]
},
{
"cell_type": "markdown",
"id": "1ab9a06f",
"metadata": {},
"source": [
"## sine-Gordon\n",
"\n",
"We now look at numerical experiments for the sine-Gordon equation: \n",
"\n",
"\\begin{align}\n",
" u_{tt} - c^2 u_{xx} + \\sin u = 0.\n",
"\\end{align}\n",
"\n",
"This is a *nonlinear equation*. "
]
},
{
"cell_type": "code",
"execution_count": 29,
"id": "8e71f4ad",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"σ = 0.5\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"\u001b[36m\u001b[1m[ \u001b[22m\u001b[39m\u001b[36m\u001b[1mInfo: \u001b[22m\u001b[39mSaved animation to c:\\Users\\math5\\Math 5485\\Math5486\\pics\\sine-Gordon.gif\n"
]
},
{
"data": {
"text/html": [
"
"
],
"text/plain": [
"Plots.AnimatedGif(\"c:\\\\Users\\\\math5\\\\Math 5485\\\\Math5486\\\\pics\\\\sine-Gordon.gif\")"
]
},
"execution_count": 29,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# | code-fold: false\n",
"\n",
"function fdm_sineGordon(m, n, c, u0, v0; a=0., b=1., g_a=0., g_b=0., T=.2)\n",
" # space discretisation\n",
" h = (b-a)/m\n",
" x = collect(a:h:b)\n",
"\n",
" # time discretisation\n",
" τ = T/n \n",
" t = collect(0:τ:T)\n",
"\n",
" σ = c*τ/h\n",
" println(\"σ = $σ\")\n",
"\n",
" u = zeros(m+1, n+1)\n",
"\n",
" # BCs \n",
" u[1,:] .= typeof(g_a) <: Function ? g_a.(t) : g_a\n",
" u[m+1,:] .= typeof(g_b) <: Function ? g_b.(t) : g_b\n",
"\n",
" # IC\n",
" u[:,1] = u0.(x)\n",
"\n",
" # u[•,2]\n",
" for i in 2:m\n",
" # Approximation of u_xx at t=0\n",
" u_xx = (u[i+1,1] - 2u[i,1] + u[i-1,1]) / h^2\n",
" u[i,2] = u[i,1] + τ*v0(x[i]) + (τ^2/2) * (c^2 * u_xx - sin(u[i,1]))\n",
" end\n",
"\n",
" for j ∈ 2:n\n",
" for i ∈ 2:m\n",
" u[i,j+1] = (σ^2) * u[i-1,j] + 2(1 - σ^2) * u[i,j] + (σ^2) * u[i+1,j] - u[i,j-1] - (τ^2) * sin(u[i,j])\n",
" end\n",
" end\n",
"\n",
" return x, t, u\n",
"end\n",
"\n",
"c = 1.\n",
"u0(x) = exp(-500(x-.5)^2)\n",
"#4 * atan(exp((x - 0.5) / sqrt(1 - 0.5^2))) \n",
"v0(x) = 0.\n",
"\n",
"x, t, u = fdm_sineGordon(100, 2000, c, u0, v0; T=10.0)\n",
"\n",
"anim = @animate for j ∈ 1:10:length(t)\n",
" plot(x, u[:, j]; \n",
" title = L\"Sine-Gordon Equation: $u_{tt} - c^2 u_{xx} + \\sin(u) = 0$\",\n",
" label = \"t = $(round(t[j], digits=2))\", \n",
" legend = :topright,\n",
" ylims = [-1, 1], \n",
" lw = 2)\n",
"end\n",
"gif(anim, \"pics/sine-Gordon.gif\", fps=15)"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "4756bc28",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"σ = 0.5\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"\u001b[36m\u001b[1m[ \u001b[22m\u001b[39m\u001b[36m\u001b[1mInfo: \u001b[22m\u001b[39mSaved animation to c:\\Users\\math5\\Math 5485\\Math5486\\pics\\sine-Gordon-2.gif\n"
]
},
{
"data": {
"text/html": [
"
"
],
"text/plain": [
"Plots.AnimatedGif(\"c:\\\\Users\\\\math5\\\\Math 5485\\\\Math5486\\\\pics\\\\sine-Gordon-2.gif\")"
]
},
"execution_count": 30,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# | code-fold: show\n",
"\n",
"c = 1.\n",
"u_stat(x) = 4 * atan(exp((x-1/2)/c))\n",
"u0(x) = u_stat(x) #+ sin(x)/10\n",
"v0(x) = 0.\n",
"\n",
"x, t, u = fdm_sineGordon(100, 100, c, u0, v0; a=-10, b=10 ,g_a=0., g_b=2π, T=10.0)\n",
"\n",
"anim = @animate for j ∈ 1:10:length(t)\n",
" plot(x, u[:, j]; \n",
" title = \"Sine-Gordon (stationary solution)\",\n",
" label = \"t = $(round(t[j], digits=2))\", \n",
" legend = :topright,\n",
" #ylims = [-1, 1], \n",
" lw = 2)\n",
" plot!(u_stat, linestyle=:dot, label=\"u_0\")\n",
"end\n",
"gif(anim, \"pics/sine-Gordon-2.gif\", fps=15)"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "c7bed05a",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"σ = 0.075\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"\u001b[36m\u001b[1m[ \u001b[22m\u001b[39m\u001b[36m\u001b[1mInfo: \u001b[22m\u001b[39mSaved animation to c:\\Users\\math5\\Math 5485\\Math5486\\pics\\sine-Gordon-3.gif\n"
]
},
{
"data": {
"text/html": [
"
"
],
"text/plain": [
"Plots.AnimatedGif(\"c:\\\\Users\\\\math5\\\\Math 5485\\\\Math5486\\\\pics\\\\sine-Gordon-3.gif\")"
]
},
"execution_count": 31,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# | code-fold: show\n",
"\n",
"c = 1.0\n",
"ω = 0.6 \n",
"η = sqrt(1 - ω^2)\n",
"\n",
"u0(x) = 0.0\n",
"v0(x) = (4 * ω * η * cosh(η * x)) / (η^2 * cosh(η * x)^2 + ω^2)\n",
"\n",
"x, t, u = fdm_sineGordon(300, 4000, c, u0, v0; \n",
" a=-20.0, b=20.0, T=40.0)\n",
"\n",
"anim = @animate for j ∈ 1:20:length(t)\n",
" plot(x, u[:, j], \n",
" title = L\"Sine-Gordon breather ($\\omega = %$ω$)\",\n",
" xlabel = \"x\", ylabel = \"u(x,t)\",\n",
" label = \"t = $(round(t[j], digits=1))\",\n",
" ylims = [-4, 4], lw = 2, color = :blue)\n",
"end\n",
"\n",
"gif(anim, \"pics/sine-Gordon-3.gif\", fps = 20)"
]
},
{
"cell_type": "markdown",
"id": "b1d5ab0c",
"metadata": {},
"source": [
":::{.callout-tip}\n",
"# TO-DO\n",
"\n",
"* A6 due a week today\n",
"* Prepare your presentations \n",
"* (Sign-up to office hours if you like)\n",
"\n",
"Previous reading:\n",
"\n",
"* Chapter 1 reading:\n",
" * @Burden sections 5.1 - 5.6, \n",
" * @FNC [zero stability](https://fncbook.com/zerostability/) \n",
"* Chapter 2 reading:\n",
" * Recap matrices: sections 7.1 - 7.2 of @Burden \n",
" * Jacobi & Gauss-Seidel: section 7.3 @Burden \n",
" * CG: section 7.6 @Burden\n",
"* Chapter 3 reading: \n",
" * Sections 11.1 - 11.2 of @Burden: shooting methods\n",
" * Lecture 16: @FNC [Galerkin method](https://fncbook.com/galerkin/)\n",
" * Lecture 17: Chapter 14.4 of @Suli: Theorem 14.6. \n",
" * Lecture 19: @Olver2014 Chapter 5\n",
" * @FNC 11.1, 11.2, 12.4 \n",
"\n",
":::"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Julia 1.11",
"language": "julia",
"name": "julia-1.11"
},
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "1.11.6"
}
},
"nbformat": 4,
"nbformat_minor": 5
}