{ "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 }