{ "cells": [ { "cell_type": "markdown", "source": [ "# Linear rotating shallow water dynamics\n", "\n", "\n", "\n", "This example solves the linear 1D rotating shallow water equations\n", "for the $u(x, t)$, $v(x, t)$ and the surface surface elevation $\\eta(x, t)$,\n", "for a fluid with constant rest-depth $H$. That is, the total fluid's depth\n", "is $H + \\eta(x, t)$ with $|\\eta| \\ll H$.\n", "\n", "The linearized equations for the evolution of $u$, $v$, $\\eta$ are:\n", "\n", "$$\n", "\\begin{aligned}\n", "\\partial_t u - f v & = - g \\partial_x \\eta - \\mathrm{D} u, \\\\\n", "\\partial_t v + f u & = - \\mathrm{D} v, \\\\\n", "\\partial_t \\eta + H \\partial_x u & = - \\mathrm{D} \\eta.\n", "\\end{aligned}\n", "$$\n", "\n", "Above, $g$ is the gravitational acceleration, $f$ is the Coriolis parameter, and\n", "$\\mathrm{D}$ indicates a hyperviscous linear operator of the form $(-1)^{n_ν} ν \\nabla^{2 n_ν}$,\n", "with $ν$ the viscosity coefficient and $n_ν$ the order of the operator.\n", "\n", "Rotation introduces the deformation length scale, $L_d = \\sqrt{g H} / f$. Disturbances with\n", "length scales much smaller than $L_d$ don't \"feel\" the rotation and propagate as inertia-gravity\n", "waves. Disturbances with length scales comparable or larger than $L_d$ should be approximately\n", "in geostrophic balance, i.e., the Coriolis acceleration $f \\widehat{\\bm{z}} \\times \\bm{u}$\n", "should be in approximate balance with the pressure gradient $-g \\bm{\\nabla} \\eta$." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "using FourierFlows, CairoMakie, Printf, Random, JLD2\n", "using LinearAlgebra: mul!, ldiv!" ], "metadata": {}, "execution_count": 1 }, { "cell_type": "markdown", "source": [ "## Coding up the equations\n", "### A demonstration of FourierFlows.jl framework\n", "\n", "What follows is a step-by-step tutorial demonstrating how you can create your own solver\n", "for an equation of your liking." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "The basic building blocks for a `FourierFlows.Problem` are:\n", "- `Grid` struct containining the physical and wavenumber grid for the problem,\n", "- `Params` struct containining all the parameters of the problem,\n", "- `Vars` struct containining arrays with the variables used in the problem,\n", "- `Equation` struct containining the coefficients of the linear operator $L$ and the function that computes the nonlinear terms, usually named `calcN!()`.\n", "\n", "The `Grid` structure is provided by FourierFlows.jl. We simply have to call one of either\n", "`OneDGrid()`, `TwoDGrid()`, or `ThreeDGrid()` constructors, depending on the dimensionality\n", "of the problem. All other structs mentioned above are problem-specific and need to be constructed\n", "for every set of equations we want to solve." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "First let's construct the `Params` struct that contains all parameters of the problem." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "struct Params{T} <: AbstractParams\n", " ν :: T # Hyperviscosity coefficient\n", " nν :: Int # Order of the hyperviscous operator\n", " g :: T # Gravitational acceleration\n", " H :: T # Fluid depth\n", " f :: T # Coriolis parameter\n", "end\n", "nothing #hide" ], "metadata": {}, "execution_count": 2 }, { "cell_type": "markdown", "source": [ "Now the `Vars` struct that contains all variables used in this problem. For this\n", "problem `Vars` includes the representations of the flow fields in physical space\n", "`u`, `v` and `η` and their Fourier transforms `uh`, `vh`, and `ηh`." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "struct Vars{Aphys, Atrans} <: AbstractVars\n", " u :: Aphys\n", " v :: Aphys\n", " η :: Aphys\n", " uh :: Atrans\n", " vh :: Atrans\n", " ηh :: Atrans\n", "end\n", "nothing #hide" ], "metadata": {}, "execution_count": 3 }, { "cell_type": "markdown", "source": [ "A constructor populates empty arrays based on the dimension of the `grid`\n", "and then creates `Vars` struct." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "\"\"\"\n", " Vars(grid)\n", "\n", "Construct the `Vars` for 1D linear shallow water dynamics based on the dimensions of the `grid` arrays.\n", "\"\"\"\n", "function Vars(grid)\n", " Dev = typeof(grid.device)\n", " T = eltype(grid)\n", "\n", " @devzeros Dev T grid.nx u v η\n", " @devzeros Dev Complex{T} grid.nkr uh vh ηh\n", "\n", " return Vars(u, v, η, uh, vh, ηh)\n", "end\n", "nothing #hide" ], "metadata": {}, "execution_count": 4 }, { "cell_type": "markdown", "source": [ "In Fourier space, the 1D linear shallow water dynamics read:\n", "\n", "$$\n", "\\begin{aligned}\n", "\\frac{\\partial \\hat{u}}{\\partial t} & = \\underbrace{ f \\hat{v} - i k g \\hat{\\eta} }_{N_u} \\; \\underbrace{- \\nu k^2 }_{L_u} \\hat{u} , \\\\\n", "\\frac{\\partial \\hat{v}}{\\partial t} & = \\underbrace{ - f \\hat{u} }_{N_v} \\; \\underbrace{- \\nu k^2 }_{L_v} \\hat{v} , \\\\\n", "\\frac{\\partial \\hat{\\eta}}{\\partial t} & = \\underbrace{ - i k H \\hat{u} }_{N_{\\eta}} \\; \\underbrace{- \\nu k^2 }_{L_{\\eta}} \\hat{\\eta} .\n", "\\end{aligned}\n", "$$\n", "Although, e.g., terms involving the Coriolis accelaration are, in principle, linear we include\n", "them in the nonlinear term $N$ because they render the linear operator $L$ non-diagonal.\n", "\n", "With these in mind, we construct function `calcN!` that computes the nonlinear terms." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "\"\"\"\n", " calcN!(N, sol, t, clock, vars, params, grid)\n", "\n", "Compute the nonlinear terms for 1D linear shallow water dynamics.\n", "\"\"\"\n", "function calcN!(N, sol, t, clock, vars, params, grid)\n", " @. vars.uh = sol[:, 1]\n", " @. vars.vh = sol[:, 2]\n", " @. vars.ηh = sol[:, 3]\n", "\n", " @. N[:, 1] = params.f * vars.vh - im * grid.kr * params.g * vars.ηh # + f v - g ∂η/∂x\n", " @. N[:, 2] = - params.f * vars.uh # - f u\n", " @. N[:, 3] = - im * grid.kr * params.H * vars.uh # - H ∂u/∂x\n", "\n", " dealias!(N, grid)\n", "\n", " return nothing\n", "end\n", "nothing #hide" ], "metadata": {}, "execution_count": 5 }, { "cell_type": "markdown", "source": [ "Next we construct the `Equation` struct:" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "\"\"\"\n", " Equation(params, grid)\n", "\n", "Construct the equation: the linear part, in this case the hyperviscous dissipation,\n", "and the nonlinear part, which is computed by `calcN!` function.\n", "\"\"\"\n", "function Equation(params, grid)\n", " T = eltype(grid)\n", " dev = grid.device\n", "\n", " L = zeros(dev, T, (grid.nkr, 3))\n", " D = @. - params.ν * grid.kr^(2*params.nν)\n", "\n", " L[:, 1] .= D # for u equation\n", " L[:, 2] .= D # for v equation\n", " L[:, 3] .= D # for η equation\n", "\n", " return FourierFlows.Equation(L, calcN!, grid)\n", "end\n", "nothing #hide" ], "metadata": {}, "execution_count": 6 }, { "cell_type": "markdown", "source": [ "We now have all necessary building blocks to construct a `FourierFlows.Problem`.\n", "It would be useful, however, to define some more \"helper functions\". For example,\n", "a function that updates all variables given the solution `sol` which comprises $\\hat{u}$,\n", "$\\hat{v}$ and $\\hat{\\eta}$:" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "\"\"\"\n", " updatevars!(prob)\n", "\n", "Update the variables in `prob.vars` using the solution in `prob.sol`.\n", "\"\"\"\n", "function updatevars!(prob)\n", " vars, grid, sol = prob.vars, prob.grid, prob.sol\n", "\n", " @. vars.uh = sol[:, 1]\n", " @. vars.vh = sol[:, 2]\n", " @. vars.ηh = sol[:, 3]\n", "\n", " ldiv!(vars.u, grid.rfftplan, deepcopy(sol[:, 1])) # use deepcopy() because irfft destroys its input\n", " ldiv!(vars.v, grid.rfftplan, deepcopy(sol[:, 2])) # use deepcopy() because irfft destroys its input\n", " ldiv!(vars.η, grid.rfftplan, deepcopy(sol[:, 3])) # use deepcopy() because irfft destroys its input\n", "\n", " return nothing\n", "end\n", "nothing #hide" ], "metadata": {}, "execution_count": 7 }, { "cell_type": "markdown", "source": [ "Another useful function is one that prescribes an initial condition to the state variable `sol`." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "\"\"\"\n", " set_uvη!(prob, u0, v0, η0)\n", "\n", "Set the state variable `prob.sol` as the Fourier transforms of `u0`, `v0`, and `η0`\n", "and update all variables in `prob.vars`.\n", "\"\"\"\n", "function set_uvη!(prob, u0, v0, η0)\n", " vars, grid, sol = prob.vars, prob.grid, prob.sol\n", "\n", " A = typeof(vars.u) # determine the type of vars.u\n", "\n", " # below, e.g., A(u0) converts u0 to the same type as vars expects\n", " # (useful when u0 is a CPU array but grid.device is GPU)\n", " mul!(vars.uh, grid.rfftplan, A(u0))\n", " mul!(vars.vh, grid.rfftplan, A(v0))\n", " mul!(vars.ηh, grid.rfftplan, A(η0))\n", "\n", " @. sol[:, 1] = vars.uh\n", " @. sol[:, 2] = vars.vh\n", " @. sol[:, 3] = vars.ηh\n", "\n", " updatevars!(prob)\n", "\n", " return nothing\n", "end\n", "nothing #hide" ], "metadata": {}, "execution_count": 8 }, { "cell_type": "markdown", "source": [ "## Let's prescibe parameter values and solve the PDE\n", "\n", "We are now ready to write up a program that sets up parameter values, constructs\n", "the problem `prob`, # time steps the solutions `prob.sol` and plots it." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "## Choosing a device: CPU or GPU" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "dev = CPU() # Device (CPU/GPU)\n", "nothing # hide" ], "metadata": {}, "execution_count": 9 }, { "cell_type": "markdown", "source": [ "## Numerical parameters and time-stepping parameters" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ " nx = 512 # grid resolution\n", "stepper = \"FilteredRK4\" # timestepper\n", " dt = 20.0 # timestep (s)\n", " nsteps = 320 # total number of time-steps\n", "nothing # hide" ], "metadata": {}, "execution_count": 10 }, { "cell_type": "markdown", "source": [ "## Physical parameters" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "Lx = 500e3 # Domain length (m)\n", "g = 9.8 # Gravitational acceleration (m s⁻²)\n", "H = 200.0 # Fluid depth (m)\n", "f = 1e-2 # Coriolis parameter (s⁻¹)\n", "ν = 100.0 # Viscosity (m² s⁻¹)\n", "nν = 1 # Viscosity order (nν = 1 means Laplacian ∇²)\n", "nothing # hide" ], "metadata": {}, "execution_count": 11 }, { "cell_type": "markdown", "source": [ "## Construct the `struct`s and you are ready to go!" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Create a `grid` and also `params`, `vars`, and the `equation` structs. Then\n", "give them all as input to the `FourierFlows.Problem()` constructor to get a\n", "problem struct, `prob`, that contains all of the above." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ " grid = OneDGrid(dev; nx, Lx)\n", " params = Params(ν, nν, g, H, f)\n", " vars = Vars(grid)\n", "equation = Equation(params, grid)\n", "\n", " prob = FourierFlows.Problem(equation, stepper, dt, grid, vars, params)\n", "nothing #hide" ], "metadata": {}, "execution_count": 12 }, { "cell_type": "markdown", "source": [ "## Setting initial conditions" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "For initial condition we take the fluid at rest ($u = v = 0$). The free surface elevation\n", "is perturbed from its rest position ($\\eta=0$); the disturbance we impose a Gaussian\n", "bump with half-width greater than the deformation radius and on top of that we\n", "superimpose some random noise with scales smaller than the deformation radius.\n", "We mask the small-scale perturbations so that it only applies in the central part\n", "of the domain by applying\n", "\n", "The system develops geostrophically-balanced jets around the Gaussian bump,\n", "while the smaller-scale noise propagates away as inertia-gravity waves." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "First let's construct the Gaussian bump." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "gaussian_width = 6e3\n", "gaussian_amplitude = 3.0\n", "gaussian_bump = @. gaussian_amplitude * exp( - grid.x^2 / (2*gaussian_width^2) )\n", "\n", "fig = Figure(size = (600, 260))\n", "ax = Axis(fig[1, 1];\n", " xlabel = \"x [km]\",\n", " ylabel = \"η [m]\",\n", " title = \"A gaussian bump with half-width ≈ \" * string(gaussian_width/1e3) * \" km\",\n", " limits = ((-Lx/2e3, Lx/2e3), nothing))\n", "\n", "lines!(ax, grid.x/1e3, gaussian_bump; # divide with 1e3 to convert m -> km\n", " color = (:black, 0.7),\n", " linewidth = 2)\n", "\n", "save(\"gaussian_bump.svg\", fig); nothing # hide" ], "metadata": {}, "execution_count": 13 }, { "cell_type": "markdown", "source": [ "![](gaussian_bump.svg)" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Next the noisy perturbation. The `mask` is simply a product of hyperbolic tangent functions." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "mask = @. 1/4 * (1 + tanh( -(grid.x - 100e3) / 10e3)) * (1 + tanh( (grid.x + 100e3) / 10e3))\n", "\n", "noise_amplitude = 0.1 # the amplitude of the noise for η(x, t=0) (m)\n", "η_noise = noise_amplitude * Random.randn(size(grid.x))\n", "@. η_noise *= mask # mask the noise\n", "\n", "fig = Figure(size = (600, 520))\n", "\n", "kwargs = (xlabel = \"x [km]\", limits = ((-Lx/2e3, Lx/2e3), nothing))\n", "\n", "ax1 = Axis(fig[1, 1]; ylabel = \"η [m]\", title = \"small-scale noise\", kwargs...)\n", "\n", "ax2 = Axis(fig[2, 1]; ylabel = \"mask\", kwargs...)\n", "\n", "lines!(ax1, grid.x/1e3, η_noise; # divide with 1e3 to convert m -> km\n", " color = (:black, 0.7),\n", " linewidth = 3)\n", "\n", "lines!(ax2, grid.x/1e3, mask; # divide with 1e3 to convert m -> km\n", " color = (:gray, 0.7),\n", " linewidth = 2)\n", "\n", "save(\"noise-mask.svg\", fig); nothing # hide" ], "metadata": {}, "execution_count": 14 }, { "cell_type": "markdown", "source": [ "![](noise-mask.svg)" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Sum the Gaussian bump and the noise and then call `set_uvη!()` to set the initial condition to the problem `prob`." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "η0 = @. gaussian_bump + η_noise\n", "u0 = zeros(grid.nx)\n", "v0 = zeros(grid.nx)\n", "\n", "set_uvη!(prob, u0, v0, η0)\n", "\n", "fig = Figure(size = (600, 260))\n", "\n", "ax = Axis(fig[1, 1];\n", " xlabel = \"x [km]\",\n", " ylabel = \"η [m]\",\n", " title = \"initial surface elevation, η(x, t=0)\",\n", " limits = ((-Lx/2e3, Lx/2e3), nothing))\n", "\n", "lines!(ax, grid.x/1e3, η0; # divide with 1e3 to convert m -> km\n", " color = (:black, 0.7),\n", " linewidth = 2)\n", "\n", "save(\"initial_eta.svg\", fig); nothing # hide" ], "metadata": {}, "execution_count": 15 }, { "cell_type": "markdown", "source": [ "![](initial_eta.svg)" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "## Saving output" ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "Output\n ├──── prob: FourierFlows.Problem{DataType, Matrix{ComplexF64}, Float64, Matrix{Float64}}\n ├──── path: ./linear_swe.jld2\n └── fields: Dict{Symbol, Function}(:sol => Main.var\"##226\".get_sol)" }, "metadata": {}, "execution_count": 16 } ], "cell_type": "code", "source": [ "filepath = \".\"\n", "filename = joinpath(filepath, \"linear_swe.jld2\")\n", "\n", "get_sol(prob) = prob.sol\n", "\n", "out = Output(prob, filename, (:sol, get_sol))" ], "metadata": {}, "execution_count": 16 }, { "cell_type": "markdown", "source": [ "We call `saveproblem` to we write the problem's configuration parameters to the .jld2 file." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "saveproblem(out)" ], "metadata": {}, "execution_count": 17 }, { "cell_type": "markdown", "source": [ "## Time-stepping the `Problem` forward" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "for j = 0:nsteps\n", " updatevars!(prob)\n", " stepforward!(prob)\n", " saveoutput(out)\n", "end" ], "metadata": {}, "execution_count": 18 }, { "cell_type": "markdown", "source": [ "## Visualizing the simulation" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "First we load the saved output files." ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "512-element Vector{Float64}:\n -250000.0\n -249023.4375\n -248046.875\n -247070.3125\n -246093.75\n -245117.1875\n -244140.625\n -243164.0625\n -242187.5\n -241210.9375\n ⋮\n 241210.9375\n 242187.5\n 243164.0625\n 244140.625\n 245117.1875\n 246093.75\n 247070.3125\n 248046.875\n 249023.4375" }, "metadata": {}, "execution_count": 19 } ], "cell_type": "code", "source": [ "using JLD2\n", "\n", "file = jldopen(out.path)\n", "iterations = parse.(Int, keys(file[\"snapshots/t\"]))\n", "\n", "nx = file[\"grid/nx\"]\n", " x = file[\"grid/x\"]" ], "metadata": {}, "execution_count": 19 }, { "cell_type": "markdown", "source": [ "Then we animate the output. We use Makie's `Observable` to animate the data. To dive into how\n", "`Observable`s work we refer to [Makie.jl's Documentation](https://makie.juliaplots.org/stable/documentation/nodes/index.html)." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ " n = Observable(1)\n", "\n", "u = @lift irfft(file[string(\"snapshots/sol/\", iterations[$n])][:, 1], nx)\n", "v = @lift irfft(file[string(\"snapshots/sol/\", iterations[$n])][:, 2], nx)\n", "η = @lift irfft(file[string(\"snapshots/sol/\", iterations[$n])][:, 3], nx)\n", "\n", "toptitle = @lift \"t = \" * @sprintf(\"%.1f\", file[string(\"snapshots/t/\", iterations[$n])]/60) * \" min\"\n", "\n", "fig = Figure(size = (600, 800))\n", "\n", "kwargs_η = (xlabel = \"x [km]\", limits = ((-Lx/2e3, Lx/2e3), nothing))\n", "kwargs_uv = (xlabel = \"x [km]\", limits = ((-Lx/2e3, Lx/2e3), (-0.3, 0.3)))\n", "\n", "ax_η = Axis(fig[2, 1]; ylabel = \"η [m]\", title = toptitle, kwargs_η...)\n", "\n", "ax_u = Axis(fig[3, 1]; ylabel = \"u [m s⁻¹]\", kwargs_uv...)\n", "\n", "ax_v = Axis(fig[4, 1]; ylabel = \"v [m s⁻¹]\", kwargs_uv...)\n", "\n", "Ld = @sprintf \"%.2f\" sqrt(g * H) / f /1e3 # divide with 1e3 to convert m -> km\n", "title = \"Deformation radius √(gh) / f = \"*string(Ld)*\" km\"\n", "\n", "fig[1, 1] = Label(fig, title, fontsize=24, tellwidth=false)\n", "\n", "lines!(ax_η, grid.x/1e3, η; # divide with 1e3 to convert m -> km\n", " color = (:blue, 0.7))\n", "\n", "lines!(ax_u, grid.x/1e3, u; # divide with 1e3 to convert m -> km\n", " color = (:red, 0.7))\n", "\n", "lines!(ax_v, grid.x/1e3, v; # divide with 1e3 to convert m -> km\n", " color = (:green, 0.7))\n", "\n", "frames = 1:length(iterations)\n", "\n", "record(fig, \"onedshallowwater.mp4\", frames, framerate=18) do i\n", " n[] = i\n", "end\n", "nothing #hide" ], "metadata": {}, "execution_count": 20 }, { "cell_type": "markdown", "source": [ "![](onedshallowwater.mp4)" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "## Geostrophic balance" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "It is instructive to compare the solution for $\\bm{u}$ with its geostrophically balanced\n", "approximation, $f \\widehat{\\bm{z}} \\times \\bm{u}_{\\rm geostrophic} = - g \\bm{\\nabla} \\eta$, i.e.,\n", "\n", "$$\n", "\\begin{aligned}\n", "v_{\\rm geostrophic} & = \\frac{g}{f} \\frac{\\partial \\eta}{\\partial x} \\ , \\\\\n", "u_{\\rm geostrophic} & = - \\frac{g}{f} \\frac{\\partial \\eta}{\\partial y} = 0 \\ .\n", "\\end{aligned}\n", "$$" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "u_geostrophic = zeros(grid.nx) # -g/f ∂η/∂y = 0\n", "v_geostrophic = params.g / params.f * irfft(im * grid.kr .* vars.ηh, grid.nx) #g/f ∂η/∂x\n", "\n", "nothing # hide" ], "metadata": {}, "execution_count": 21 }, { "cell_type": "markdown", "source": [ "The geostrophic solution should capture well the the behavior of the flow in the center\n", "of the domain, after small-scale disturbances propagate away. Let's plot and see!" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "fig = Figure(size = (600, 600))\n", "\n", "kwargs = (xlabel = \"x [km]\", limits = ((-Lx/2e3, Lx/2e3), (-0.3, 0.3)))\n", "\n", "ax_u = Axis(fig[2, 1]; ylabel = \"u [m s⁻¹]\", kwargs...)\n", "\n", "ax_v = Axis(fig[3, 1]; ylabel = \"v [m s⁻¹]\", kwargs...)\n", "\n", "fig[1, 1] = Label(fig, \"Geostrophic balance\", fontsize=24, tellwidth=false)\n", "\n", "lines!(ax_u, grid.x/1e3, vars.u; # divide with 1e3 to convert m -> km\n", " label = \"u\",\n", " linewidth = 3,\n", " color = (:red, 0.7))\n", "\n", "lines!(ax_u, grid.x/1e3, u_geostrophic; # divide with 1e3 to convert m -> km\n", " label = \"- g/f ∂η/∂y\",\n", " linewidth = 3,\n", " color = (:purple, 0.7))\n", "\n", "axislegend(ax_u)\n", "\n", "lines!(ax_v, grid.x/1e3, vars.v; # divide with 1e3 to convert m -> km\n", " label = \"v\",\n", " linewidth = 3,\n", " color = (:green, 0.7))\n", "\n", "lines!(ax_v, grid.x/1e3, v_geostrophic; # divide with 1e3 to convert m -> km\n", " label = \"g/f ∂η/∂x\",\n", " linewidth = 3,\n", " color = (:purple, 0.7))\n", "\n", "axislegend(ax_v)\n", "\n", "save(\"geostrophic_balance.svg\", fig); nothing # hide" ], "metadata": {}, "execution_count": 22 }, { "cell_type": "markdown", "source": [ "![](geostrophic_balance.svg)" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "---\n", "\n", "*This notebook was generated using [Literate.jl](https://github.com/fredrikekre/Literate.jl).*" ], "metadata": {} } ], "nbformat_minor": 3, "metadata": { "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.10.0" }, "kernelspec": { "name": "julia-1.10", "display_name": "Julia 1.10.0", "language": "julia" } }, "nbformat": 4 }