From 467e82f0592309edaa687f58c9fe8e9f167d2290 Mon Sep 17 00:00:00 2001 From: "Navid C. Constantinou" Date: Sat, 27 Jul 2024 13:04:15 +1000 Subject: [PATCH 1/7] Update papers using GeophysicalFlows (#364) * Update papers using GeophysicalFlows * cleanup --- docs/src/index.md | 2 ++ 1 file changed, 2 insertions(+) diff --git a/docs/src/index.md b/docs/src/index.md index fb086f33..6abfdb27 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -84,6 +84,8 @@ The bibtex entry for the paper is: ## Papers using `GeophysicalFlows.jl` +1. Lobo, M., Griffies, S. M., and Zhang, W. (2024) Vertical structure of baroclinic instability in a three-layer quasi-geostrophic model over a sloping bottom. _ESS Open Archive_, doi:[10.22541/essoar.172166435.57577370/v1](https://doi.org/10.22541/essoar.172166435.57577370/v1). + 1. Pudig, M. and Smith, K. S. (2024) Baroclinic turbulence above rough topography: The vortex gas and topographic turbulence regimes. _ESS Open Archive_, doi:[10.22541/essoar.171995116.60993353/v1](https://doi.org/10.22541/essoar.171995116.60993353/v1). 1. Shokar, I. J. S., Haynes, P. H. and Kerswell, R. R. (2024) Extending deep learning emulation across parameter regimes to assess stochastically driven spontaneous transition events. In ICLR 2024 Workshop on AI4DifferentialEquations in Science. url: [https://openreview.net/forum?id=7a5gUX4e5q](https://openreview.net/forum?id=7a5gUX4e5q). From 14f581bd6d36e6ba74cf2136eaaa1849d0c5f045 Mon Sep 17 00:00:00 2001 From: mncrowe <87329513+mncrowe@users.noreply.github.com> Date: Mon, 29 Jul 2024 15:19:30 +0100 Subject: [PATCH 2/7] Add background zonal flow `U` or `U(y)` in the `SingleLayerQG` module (#360) MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit * Update singlelayerqg.jl to include background flow * Fixed typo in previous change * added background flow test * corrected typo * Added some spaces to cheer up Greg Co-authored-by: Gregory L. Wagner * more spaces * added variable U functionality * reduced memory of Uy assignment and updated function descriptions * update function descriptions * T(U) -> convert(T, U) Co-authored-by: Gregory L. Wagner * clarify U array formulation Co-authored-by: Gregory L. Wagner * fix bug with previous push * add additional test function for U as an Array and include tests in * simplify background U and use Rossby wave test for U that does not vary in y * homogenize notation/docs * bump patch release * fix test for gpu * == -> === * use different calcN_advection! methods for U const and y-varying * minor tweaks * update singlelayergq module in docs * homogenize notation/wording * add empty line * code alignment * add Qx and Qy in params; refactor calcN_advection * clarify when U=const * some fixes + cleanup * add nonlinear advection test for finite deformation + topo * fixes nonlinear advection w U(y) + add tests * add comment for β term * remove debug statements + fix test for julia 1.6 * more tests for the nonlinear terms * fix test_1layerqg_nonlinearadvection on gpu * fix test_1layerqg_nonlinearadvection on gpu * improve docs * drop GeophysicalFlows from docs/Project.toml --------- Co-authored-by: Gregory L. Wagner Co-authored-by: Navid C. Constantinou --- Project.toml | 2 +- docs/src/modules/multilayerqg.md | 6 +- docs/src/modules/singlelayerqg.md | 37 ++++-- src/multilayerqg.jl | 5 +- src/singlelayerqg.jl | 177 +++++++++++++++++++++------- test/runtests.jl | 37 +++--- test/test_singlelayerqg.jl | 188 +++++++++++++++++++++++------- 7 files changed, 340 insertions(+), 112 deletions(-) diff --git a/Project.toml b/Project.toml index a49c334c..452978c4 100644 --- a/Project.toml +++ b/Project.toml @@ -2,7 +2,7 @@ name = "GeophysicalFlows" uuid = "44ee3b1c-bc02-53fa-8355-8e347616e15e" license = "MIT" authors = ["Navid C. Constantinou ", "Gregory L. Wagner ", "and contributors"] -version = "0.16.1" +version = "0.16.2" [deps] CUDA = "052768ef-5323-5732-b1bb-66c8b64840ba" diff --git a/docs/src/modules/multilayerqg.md b/docs/src/modules/multilayerqg.md index 96d0488b..d9717c8c 100644 --- a/docs/src/modules/multilayerqg.md +++ b/docs/src/modules/multilayerqg.md @@ -71,10 +71,14 @@ with ```math \begin{aligned} \partial_y Q_j &\equiv \beta - \partial_y^2 U_j - (1-\delta_{j,1}) F_{j-1/2, j} (U_{j-1} - U_j) - (1 - \delta_{j,n}) F_{j+1/2, j} (U_{j+1} - U_j) + \delta_{j,n} \partial_y \eta , \\ -\partial_x Q_j &\equiv \delta_{j, n} \partial_x \eta . +\partial_x Q_j & \equiv \delta_{j, n} \partial_x \eta , \end{aligned} ``` +the background PV gradient components in each layer and with +``\mathsf{J}(a, b) = (\partial_x a)(\partial_y b)-(\partial_y a)(\partial_x b)`` is the +two-dimensional Jacobian. On the right hand side, ``\mu`` is linear bottom drag, and ``\nu`` is +hyperviscosity of order ``n_\nu``. Plain old viscosity corresponds to ``n_\nu = 1``. ### Implementation diff --git a/docs/src/modules/singlelayerqg.md b/docs/src/modules/singlelayerqg.md index 6f6e382c..e8e09603 100644 --- a/docs/src/modules/singlelayerqg.md +++ b/docs/src/modules/singlelayerqg.md @@ -20,33 +20,50 @@ radius follows is said to obey equivalent-barotropic dynamics. We denote the sum vorticity and the vortex stretching contributions to the QGPV with ``q \equiv \nabla^2 \psi - \psi / \ell^2``. Also, we denote the topographic PV with ``\eta \equiv f_0 h / H``. -The dynamical variable is ``q``. Thus, the equation solved by the module is: +The dynamical variable is ``q``. Including an imposed zonal flow ``U(y)``, the equation of motion is: ```math -\partial_t q + \mathsf{J}(\psi, q + \eta) + \beta \partial_x \psi = -\underbrace{-\left[\mu + \nu(-1)^{n_\nu} \nabla^{2n_\nu} \right] q}_{\textrm{dissipation}} + F . +\partial_t q + \mathsf{J}(\psi, q) + (U - \partial_y\psi) \partial_x Q + U \partial_x q + (\partial_y Q)(\partial_x \psi) = \underbrace{-\left[\mu + \nu(-1)^{n_\nu} \nabla^{2n_\nu} \right] q}_{\textrm{dissipation}} + F , ``` -where ``\mathsf{J}(a, b) = (\partial_x a)(\partial_y b)-(\partial_y a)(\partial_x b)`` is the -two-dimensional Jacobian. On the right hand side, ``F(x, y, t)`` is forcing, ``\mu`` is +with + +```math +\begin{aligned} +\partial_y Q &\equiv \beta - \partial_y^2 U + \partial_y \eta , \\ +\partial_x Q &\equiv \partial_x \eta , +\end{aligned} +``` + +the background PV gradient components, and with +``\mathsf{J}(a, b) = (\partial_x a)(\partial_y b) - (\partial_y a)(\partial_x b)`` +the two-dimensional Jacobian. On the right hand side, ``F(x, y, t)`` is forcing, ``\mu`` is linear drag, and ``\nu`` is hyperviscosity of order ``n_\nu``. Plain old viscosity corresponds to ``n_\nu = 1``. +In the case that the imposed background zonal flow is just a constant, the above simplifies to: + +```math +\partial_t q + \mathsf{J}(\psi, q + \eta) + U \partial_x (q + \eta) + β \partial_x \psi = \underbrace{-\left[\mu + \nu(-1)^{n_\nu} \nabla^{2n_\nu} \right] q}_{\textrm{dissipation}} + F , +``` + +and thus the advection of ``q`` can be incorporated in the linear term ``L``. ### Implementation The equation is time-stepped forward in Fourier space: ```math -\partial_t \widehat{q} = - \widehat{\mathsf{J}(\psi, q + \eta)} + \beta \frac{i k_x}{|𝐤|^2 + 1/\ell^2} \widehat{q} - \left(\mu + \nu |𝐤|^{2n_\nu} \right) \widehat{q} + \widehat{F} . +\partial_t \widehat{q} = - \widehat{\mathsf{J}(\psi, q)} - \widehat{U \partial_x Q} - \widehat{U \partial_x q} ++ \widehat{(\partial_y \psi) (\partial_x Q)} - \widehat{(\partial_x \psi)(\partial_y Q)} - \left(\mu + \nu |𝐤|^{2n_\nu} \right) \widehat{q} + \widehat{F} . ``` +In doing so the Jacobian is computed in the conservative form: ``\mathsf{J}(f,g) = +\partial_y [ (\partial_x f) g] - \partial_x[ (\partial_y f) g]``. + The state variable `sol` is the Fourier transform of the sum of relative vorticity and vortex stretching (when the latter is applicable), [`qh`](@ref GeophysicalFlows.SingleLayerQG.Vars). -The Jacobian is computed in the conservative form: ``\mathsf{J}(f, g) = -\partial_y [ (\partial_x f) g] - \partial_x[ (\partial_y f) g]``. - The linear operator is constructed in `Equation` ```@docs @@ -122,4 +139,4 @@ Other diagnostic include: [`energy_dissipation`](@ref GeophysicalFlows.SingleLay (barotropic quasi-geostrophic flow with ``\beta=0``) above topography. - [`examples/singlelayerqg_decaying_barotropic_equivalentbarotropic.jl`](@ref singlelayerqg_decaying_barotropic_equivalentbarotropic_example): - Simulate two dimensional turbulence (``\beta=0``) with both infinite and finite Rossby radius of deformation and compares the evolution of the two. \ No newline at end of file + Simulate two dimensional turbulence (``\beta=0``) with both infinite and finite Rossby radius of deformation and compares the evolution of the two. diff --git a/src/multilayerqg.jl b/src/multilayerqg.jl index 4e35d679..99668a8b 100644 --- a/src/multilayerqg.jl +++ b/src/multilayerqg.jl @@ -66,7 +66,7 @@ Keyword arguments - `Ly`: Extent of the ``y``-domain. - `f₀`: Constant planetary vorticity. - `β`: Planetary vorticity ``y``-gradient. - - `U`: The imposed constant zonal flow ``U(y)`` in each fluid layer. + - `U`: Imposed background constant zonal flow ``U(y)`` in each fluid layer. - `H`: Rest height of each fluid layer. - `b`: Boussinesq buoyancy of each fluid layer. - `eta`: Periodic component of the topographic potential vorticity. @@ -757,7 +757,7 @@ function calcN_advection!(N, sol, vars, params, grid) @. vars.u += params.U # add the imposed zonal flow U uQx, uQxh = vars.q, vars.uh # use vars.q and vars.uh as scratch variables - @. uQx = vars.u * params.Qx # (U+u)*∂Q/∂x + @. uQx = vars.u * params.Qx # (U+u)*∂Q/∂x fwdtransform!(uQxh, uQx, params) @. N = - uQxh # -\hat{(U+u)*∂Q/∂x} @@ -803,6 +803,7 @@ function calcN_linearadvection!(N, sol, vars, params, grid) @. vars.vh = im * grid.kr * vars.ψh invtransform!(vars.u, vars.uh, params) + @. vars.u += params.U # add the imposed zonal flow U uQx, uQxh = vars.q, vars.uh # use vars.q and vars.uh as scratch variables @. uQx = vars.u * params.Qx # (U+u)*∂Q/∂x diff --git a/src/singlelayerqg.jl b/src/singlelayerqg.jl index 2ce0fc87..baaf26e0 100644 --- a/src/singlelayerqg.jl +++ b/src/singlelayerqg.jl @@ -1,5 +1,7 @@ module SingleLayerQG +using LinearAlgebra + export Problem, streamfunctionfrompv!, @@ -37,6 +39,7 @@ nothingfunction(args...) = nothing Ly = Lx, β = 0.0, deformation_radius = Inf, + U = 0.0, eta = nothing, ν = 0.0, nν = 1, @@ -62,6 +65,7 @@ Keyword arguments - `Ly`: Extent of the ``y``-domain. - `β`: Planetary vorticity ``y``-gradient. - `deformation_radius`: Rossby radius of deformation; set `Inf` for purely barotropic. + - `U`: Imposed background constant zonal flow ``U(y)``. - `eta`: Topographic potential vorticity. - `ν`: Small-scale (hyper)-viscosity coefficient. - `nν`: (Hyper)-viscosity order, `nν```≥ 1``. @@ -82,6 +86,7 @@ function Problem(dev::Device=CPU(); # Physical parameters β = 0.0, deformation_radius = Inf, + U = 0.0, eta = nothing, # Drag and (hyper-)viscosity ν = 0.0, @@ -98,12 +103,15 @@ function Problem(dev::Device=CPU(); # the grid grid = TwoDGrid(dev; nx, Lx, ny, Ly, aliased_fraction, T) - x, y = gridpoints(grid) # topographic PV - eta === nothing && (eta = zeros(dev, T, (nx, ny))) + eta isa Nothing && (eta = zeros(dev, T, (nx, ny))) + + U = U isa Number ? convert(T, U) : U - params = deformation_radius == Inf ? BarotropicQGParams(grid, T(β), eta, T(μ), T(ν), nν, calcF) : EquivalentBarotropicQGParams(grid, T(β), T(deformation_radius), eta, T(μ), T(ν), nν, calcF) + params = (deformation_radius == Inf || + deformation_radius === nothing) ? BarotropicQGParams(grid, T(β), U, eta, T(μ), T(ν), nν, calcF) : + EquivalentBarotropicQGParams(grid, T(deformation_radius), T(β), U, eta, T(μ), T(ν), nν, calcF) vars = calcF == nothingfunction ? DecayingVars(grid) : (stochastic ? StochasticForcedVars(grid) : ForcedVars(grid)) @@ -120,17 +128,19 @@ end abstract type SingleLayerQGParams <: AbstractParams end """ - struct Params{T, Aphys, Atrans, ℓ} <: SingleLayerQGParams + struct Params{T, Aphys, Atrans, Tℓ, TU} <: SingleLayerQGParams The parameters for the `SingleLayerQG` problem. $(TYPEDFIELDS) """ -struct Params{T, Aphys, Atrans, ℓ} <: SingleLayerQGParams +struct Params{T, Aphys, Atrans, Tℓ, TU} <: SingleLayerQGParams "planetary vorticity ``y``-gradient" β :: T "Rossby radius of deformation" - deformation_radius :: ℓ + deformation_radius :: Tℓ + "Background flow in ``x`` direction" + U :: TU "topographic potential vorticity" eta :: Aphys "Fourier transform of topographic potential vorticity" @@ -143,31 +153,59 @@ struct Params{T, Aphys, Atrans, ℓ} <: SingleLayerQGParams nν :: Int "function that calculates the Fourier transform of the forcing, ``F̂``" calcF! :: Function + # derived params + "array containing ``x``-gradient of PV due to eta" + Qx :: Aphys + "array containing ``y``-gradient of PV due to ``U`` and topographic PV" + Qy :: Aphys end -const BarotropicQGParams = Params{<:AbstractFloat, <:AbstractArray, <:AbstractArray, Nothing} -const EquivalentBarotropicQGParams = Params{<:AbstractFloat, <:AbstractArray, <:AbstractArray, <:AbstractFloat} +const BarotropicQGParams = Params{<:AbstractFloat, <:AbstractArray, <:AbstractArray, <:Nothing, <:Any} +const EquivalentBarotropicQGParams = Params{<:AbstractFloat, <:AbstractArray, <:AbstractArray, <:AbstractFloat, <:Any} + +const SingleLayerQGconstantUParams = Params{<:AbstractFloat, <:AbstractArray, <:AbstractArray, <:Any, <:Number} +const SingleLayerQGvaryingUParams = Params{<:AbstractFloat, <:AbstractArray, <:AbstractArray, <:Any, <:AbstractArray} """ - EquivalentBarotropicQGParams(grid, β, deformation_radius, eta, μ, ν, nν, calcF) + EquivalentBarotropicQGParams(grid, deformation_radius, β, U, eta, μ, ν, nν, calcF) Return the parameters for an Equivalent Barotropic QG problem (i.e., with finite Rossby radius of deformation). """ -function EquivalentBarotropicQGParams(grid::AbstractGrid{T, A}, β, deformation_radius, eta, μ, ν, nν::Int, calcF) where {T, A} - eta_on_grid = typeof(eta) <: AbstractArray ? A(eta) : FourierFlows.on_grid(eta, grid) - etah_on_grid = rfft(eta_on_grid) +function EquivalentBarotropicQGParams(grid::AbstractGrid{T, A}, deformation_radius, β, U, eta, μ, ν, nν::Int, calcF) where {T, A} + + if U isa AbstractArray && length(U) == grid.ny + U = repeat(reshape(U, (1, grid.ny)), outer=(grid.nx, 1)) # convert to 2D + U = A(U) + end + + eta_on_grid = eta isa AbstractArray ? A(eta) : FourierFlows.on_grid(eta, grid) + etah = rfft(eta_on_grid) + + Qx = irfft(im * grid.kr .* etah, grid.nx) # ∂η/∂x + Qy = irfft(im * grid.l .* etah, grid.nx) # ∂η/∂y + + if U isa AbstractArray + Uh = rfft(U) + Uyy = irfft(- grid.l.^2 .* Uh, grid.nx) # ∂²U/∂y² + Qy .-= Uyy # -∂²U/∂y² + end + + # Note: The β-term in Qy is included in the linear term L of Equation. - return Params(β, deformation_radius, eta_on_grid, etah_on_grid, μ, ν, nν, calcF) + return Params(β, deformation_radius, U, eta_on_grid, etah, μ, ν, nν, calcF, Qx, Qy) end """ - BarotropicQGParams(grid, β, eta, μ, ν, nν, calcF) + BarotropicQGParams(grid, β, U, eta, μ, ν, nν, calcF) Return the parameters for a Barotropic QG problem (i.e., with infinite Rossby radius of deformation). """ -BarotropicQGParams(grid::AbstractGrid, β, eta, μ, ν, nν::Int, calcF) = - EquivalentBarotropicQGParams(grid, β, nothing, eta, μ, ν, nν, calcF) - +function BarotropicQGParams(grid, β, U, eta, μ, ν, nν::Int, calcF) + deformation_radius = nothing + + return EquivalentBarotropicQGParams(grid, deformation_radius, β, U, eta, μ, ν, nν, calcF) +end + # --------- # Equations @@ -175,41 +213,41 @@ BarotropicQGParams(grid::AbstractGrid, β, eta, μ, ν, nν::Int, calcF) = """ Equation(params::BarotropicQGParams, grid) + Equation(params::EquivalentBarotropicQGParams, grid) -Return the equation for a barotropic QG problem with `params` and `grid`. Linear operator -``L`` includes bottom drag ``μ``, (hyper)-viscosity of order ``n_ν`` with coefficient ``ν`` -and the ``β`` term: +Return the equation for a `SingleLayerQG` problem with `params` and `grid`. +Linear operator ``L`` includes bottom drag ``μ``, (hyper)-viscosity of order ``n_ν`` with +coefficient ``ν``, and the ``β`` term. If there is a constant background flow ``U`` that +does not vary in ``y`` then the linear term ``L`` includes also the mean advection term +by ``U``, namely ``-i k_x U```. That is: ```math -L = - μ - ν |𝐤|^{2 n_ν} + i β k_x / |𝐤|² . +L = -μ - ν |𝐤|^{2 n_ν} + i β k_x / (|𝐤|² + 1/ℓ²) - i k_x U . ``` The nonlinear term is computed via `calcN!` function. """ -function Equation(params::BarotropicQGParams, grid::AbstractGrid) +function Equation(params::BarotropicQGParams, grid) L = @. - params.μ - params.ν * grid.Krsq^params.nν + im * params.β * grid.kr * grid.invKrsq + + if params.U isa Number + @. L -= im * params.U * grid.kr + end + CUDA.@allowscalar L[1, 1] = 0 - + return FourierFlows.Equation(L, calcN!, grid) end -""" - Equation(params::EquivalentBarotropicQGParams, grid) - -Return the equation for an equivalent-barotropic QG problem with `params` and `grid`. -Linear operator ``L`` includes bottom drag ``μ``, (hyper)-viscosity of order ``n_ν`` with -coefficient ``ν`` and the ``β`` term: +function Equation(params::EquivalentBarotropicQGParams, grid) + L = @. - params.μ - params.ν * grid.Krsq^params.nν + im * params.β * grid.kr / (grid.Krsq + 1 / params.deformation_radius^2) -```math -L = -μ - ν |𝐤|^{2 n_ν} + i β k_x / (|𝐤|² + 1/ℓ²) . -``` + if params.U isa Number + @. L -= im * params.U * grid.kr + end -The nonlinear term is computed via `calcN!` function. -""" -function Equation(params::EquivalentBarotropicQGParams, grid::AbstractGrid) - L = @. - params.μ - params.ν * grid.Krsq^params.nν + im * params.β * grid.kr / (grid.Krsq + 1 / params.deformation_radius^2) CUDA.@allowscalar L[1, 1] = 0 - + return FourierFlows.Equation(L, calcN!, grid) end @@ -304,16 +342,19 @@ end # ------- """ - calcN_advection!(N, sol, t, clock, vars, params, grid) + calcN_advection!(N, sol, t, clock, vars, params::SingleLayerQGconstantUParams, grid) -Calculate the Fourier transform of the advection term, ``- 𝖩(ψ, q+η)`` in conservative -form, i.e., ``- ∂_x[(∂_y ψ)(q+η)] - ∂_y[(∂_x ψ)(q+η)]`` and store it in `N`: +Compute the advection term and stores it in `N`. The imposed zonal flow ``U`` is either +zero or constant, in which case is incorporated in the linear terms of the equation. +Thus, the nonlinear terms is ``- 𝖩(ψ, q + η)`` in conservative form, i.e., +``- ∂_x[(∂_y ψ)(q + η)] - ∂_y[(∂_x ψ)(q + η)]``: ```math N = - \\widehat{𝖩(ψ, q + η)} = - i k_x \\widehat{u (q + η)} - i k_y \\widehat{v (q + η)} . ``` """ -function calcN_advection!(N, sol, t, clock, vars, params, grid) +function calcN_advection!(N, sol, t, clock, vars, params::SingleLayerQGconstantUParams, grid) + @. vars.qh = sol streamfunctionfrompv!(vars.ψh, vars.qh, params, grid) @. vars.uh = -im * grid.l * vars.ψh @@ -335,6 +376,58 @@ function calcN_advection!(N, sol, t, clock, vars, params, grid) @. N = -im * grid.kr * uq_plus_ηh - im * grid.l * vq_plus_ηh # - ∂[u*(q+η)]/∂x - ∂[v*(q+η)]/∂y + @. N -= im * grid.kr * params.U * params.etah # - \hat{∂(Uη)/∂x} + + return nothing +end + +""" + calcN_advection!(N, sol, t, clock, vars, params::SingleLayerQGvaryingUParams, grid) + +Compute the advection term and stores it in `N`. The imposed zonal flow ``U(y)`` varies +with ``y`` and therefore is not taken care by the linear term `L` but rather is +incorporated in the nonlinear term `N`. + +```math +N = - \\widehat{𝖩(ψ, q + η)} - \\widehat{U ∂_x (q + η)} + \\widehat{(∂_x ψ)(∂_y² U)} . +``` +""" +function calcN_advection!(N, sol, t, clock, vars, params::SingleLayerQGvaryingUParams, grid) + + @. vars.qh = sol + + streamfunctionfrompv!(vars.ψh, vars.qh, params, grid) + + @. vars.uh = -im * grid.l * vars.ψh + @. vars.vh = im * grid.kr * vars.ψh + + ldiv!(vars.u, grid.rfftplan, vars.uh) + @. vars.u += params.U # add the imposed zonal flow U + + uQx, uQxh = vars.q, vars.uh # use vars.q and vars.uh as scratch variables + @. uQx = vars.u * params.Qx # (U+u)*∂Q/∂x + mul!(uQxh, grid.rfftplan, uQx) + @. N = - uQxh # -\hat{(U+u)*∂Q/∂x} + + ldiv!(vars.v, grid.rfftplan, vars.vh) + + vQy, vQyh = vars.q, vars.vh # use vars.q and vars.vh as scratch variables + @. vQy = vars.v * params.Qy # v*∂Q/∂y + mul!(vQyh, grid.rfftplan, vQy) + @. N -= vQyh # -\hat{v*∂Q/∂y} + + ldiv!(vars.q, grid.rfftplan, vars.qh) + + uq , vq = vars.u , vars.v # use vars.u and vars.v as scratch variables + uqh, vqh = vars.uh, vars.vh # use vars.uh and vars.vh as scratch variables + @. uq *= vars.q # (U+u)*q + @. vq *= vars.q # v*q + + mul!(uqh, grid.rfftplan, uq) + mul!(vqh, grid.rfftplan, vq) + + @. N -= im * grid.kr * uqh + im * grid.l * vqh # -\hat{∂[(U+u)q]/∂x} - \hat{∂[vq]/∂y} + return nothing end @@ -344,7 +437,7 @@ end Calculate the nonlinear term, that is the advection term and the forcing, ```math -N = - \\widehat{𝖩(ψ, q + η)} + F̂ . +N = - \\widehat{𝖩(ψ, q + η)} - \\widehat{U ∂_x (q + η)} + \\widehat{(∂_x ψ)(∂_y² U)} + F̂ . ``` """ function calcN!(N, sol, t, clock, vars, params, grid) diff --git a/test/runtests.jl b/test/runtests.jl index 6d542b55..95e193c1 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -48,32 +48,37 @@ for dev in devices @test test_twodnavierstokes_stochasticforcing_enstrophybudget(dev) @test test_twodnavierstokes_energyenstrophypalinstrophy(dev) @test test_twodnavierstokes_problemtype(dev, Float32) - @test TwoDNavierStokes.nothingfunction() == nothing + @test TwoDNavierStokes.nothingfunction() === nothing end @testset "SingleLayerQG" begin include("test_singlelayerqg.jl") - for deformation_radius in [Inf, 1.23] - @test test_1layerqg_rossbywave("ETDRK4", 1e-2, 20, dev, deformation_radius=deformation_radius) - @test test_1layerqg_rossbywave("FilteredETDRK4", 1e-2, 20, dev, deformation_radius=deformation_radius) - @test test_1layerqg_rossbywave("RK4", 1e-2, 20, dev, deformation_radius=deformation_radius) - @test test_1layerqg_rossbywave("FilteredRK4", 1e-2, 20, dev, deformation_radius=deformation_radius) - @test test_1layerqg_rossbywave("AB3", 1e-3, 200, dev, deformation_radius=deformation_radius) - @test test_1layerqg_rossbywave("FilteredAB3", 1e-3, 200, dev, deformation_radius=deformation_radius) - @test test_1layerqg_rossbywave("ForwardEuler", 1e-4, 2000, dev, deformation_radius=deformation_radius) - @test test_1layerqg_rossbywave("FilteredForwardEuler", 1e-4, 2000, dev, deformation_radius=deformation_radius) - @test test_1layerqg_problemtype(dev, Float32, deformation_radius=deformation_radius) + for deformation_radius in (Inf, 1.23), U₀ in (0, 0.3) + for (timestepper, dt, nsteps) in zip(("ETDRK4", "FilteredETDRK4", "RK4", "FilteredRK4", "AB3", "FilteredAB3", "ForwardEuler", "FilteredForwardEuler",), + (1e-2, 1e-2, 1e-2, 1e-2, 1e-3, 1e-3, 1e-4, 1e-4,), + (20, 20, 20, 20, 200, 200, 2000, 2000,)) + + nx = 64 + @test test_1layerqg_rossbywave(timestepper, dt, nsteps, dev, nx; deformation_radius, U=U₀) + @test test_1layerqg_rossbywave(timestepper, dt, nsteps, dev, nx; deformation_radius, U=U₀*ones((nx,))) + end + @test test_1layerqg_problemtype(dev, Float32; deformation_radius, U=U₀) end - @test test_1layerqg_advection(0.0005, "ForwardEuler", dev) + @test test_1layerqg_nonlinearadvection(0.0005, "ForwardEuler", dev, add_background_flow = false, add_topography = false) + @test test_1layerqg_nonlinearadvection(0.0005, "ForwardEuler", dev, add_background_flow = false, add_topography = true) + @test test_1layerqg_nonlinearadvection(0.0005, "ForwardEuler", dev, add_background_flow = true, background_flow_vary_in_y = false) + @test test_1layerqg_nonlinearadvection(0.0005, "ForwardEuler", dev, add_background_flow = true, background_flow_vary_in_y = false, add_topography = true) + @test test_1layerqg_nonlinearadvection(0.0005, "ForwardEuler", dev, add_background_flow = true, background_flow_vary_in_y = true) + @test test_1layerqg_nonlinearadvection_deformation(0.0005, "ForwardEuler", dev) @test test_streamfunctionfrompv(dev; deformation_radius=1.23) - @test test_1layerqg_energies_EquivalentBarotropicQG(dev; deformation_radius=1.23) @test test_1layerqg_energyenstrophy_BarotropicQG(dev) + @test test_1layerqg_energies_EquivalentBarotropicQG(dev; deformation_radius=1.23) @test test_1layerqg_deterministicforcing_energybudget(dev) @test test_1layerqg_stochasticforcing_energybudget(dev) @test test_1layerqg_deterministicforcing_enstrophybudget(dev) @test test_1layerqg_stochasticforcing_enstrophybudget(dev) - @test SingleLayerQG.nothingfunction() == nothing + @test SingleLayerQG.nothingfunction() === nothing @test_throws ErrorException("not implemented for finite deformation radius") test_1layerqg_energy_dissipation(dev; deformation_radius=2.23) @test_throws ErrorException("not implemented for finite deformation radius") test_1layerqg_enstrophy_dissipation(dev; deformation_radius=2.23) @test_throws ErrorException("not implemented for finite deformation radius") test_1layerqg_energy_work(dev; deformation_radius=2.23) @@ -98,7 +103,7 @@ for dev in devices @test test_bqgql_advection(0.0005, "ForwardEuler", dev) @test test_bqgql_energyenstrophy(dev) @test test_bqgql_problemtype(dev, Float32) - @test BarotropicQGQL.nothingfunction() == nothing + @test BarotropicQGQL.nothingfunction() === nothing end @testset "SurfaceQG" begin @@ -112,7 +117,7 @@ for dev in devices @test test_sqg_problemtype(dev, Float32) @test test_sqg_paramsconstructor(dev) @test test_sqg_noforcing(dev) - @test SurfaceQG.nothingfunction() == nothing + @test SurfaceQG.nothingfunction() === nothing end @testset "MultiLayerQG" begin diff --git a/test/test_singlelayerqg.jl b/test/test_singlelayerqg.jl index 2cca28e3..74289e03 100644 --- a/test/test_singlelayerqg.jl +++ b/test/test_singlelayerqg.jl @@ -1,10 +1,12 @@ +using LinearAlgebra + """ test_1layerqg_rossbywave(stepper, dt, nsteps, dev; kwargs...) Evolves a Rossby wave and compares with the analytic solution. """ -function test_1layerqg_rossbywave(stepper, dt, nsteps, dev::Device=CPU(); deformation_radius=Inf) - nx = 64 +function test_1layerqg_rossbywave(stepper, dt, nsteps, dev::Device=CPU(), nx=64; deformation_radius=Inf, U = 0) + Lx = 2π β = 2.0 μ = 0.0 @@ -19,15 +21,15 @@ function test_1layerqg_rossbywave(stepper, dt, nsteps, dev::Device=CPU(); deform eta(x, y) = 0 * x end - prob = SingleLayerQG.Problem(dev; nx=nx, Lx=Lx, eta=eta, deformation_radius=deformation_radius, β=β, μ=μ, ν=ν, stepper=stepper, dt=dt) + prob = SingleLayerQG.Problem(dev; nx=nx, Lx=Lx, eta=eta, deformation_radius=deformation_radius, β=β, U=U, μ=μ, ν=ν, stepper=stepper, dt=dt) sol, clock, vars, params, grid = prob.sol, prob.clock, prob.vars, prob.params, prob.grid x, y = gridpoints(grid) # the Rossby wave initial condition ampl = 1e-2 - kwave = 3 * 2π/grid.Lx - lwave = 2 * 2π/grid.Ly + kwave = 3 * 2π / grid.Lx + lwave = 2 * 2π / grid.Ly ω = -params.β * kwave / (kwave^2 + lwave^2 + 1 / deformation_radius^2) q₀ = @. ampl * cos(kwave * x) * cos(lwave * y) q₀h = rfft(q₀) @@ -38,7 +40,7 @@ function test_1layerqg_rossbywave(stepper, dt, nsteps, dev::Device=CPU(); deform dealias!(sol, grid) SingleLayerQG.updatevars!(prob) - q_theory = @. ampl * cos(kwave * (x - ω / kwave * clock.t)) * cos(lwave * y) + q_theory = @CUDA.allowscalar @. ampl * cos(kwave * (x - (U[1] + ω / kwave) * clock.t)) * cos(lwave * y) return isapprox(q_theory, vars.q, rtol=grid.nx * grid.ny * nsteps * 1e-12) end @@ -255,17 +257,26 @@ function test_1layerqg_deterministicforcing_enstrophybudget(dev::Device=CPU(); n end """ - test_1layerqg_nonlinearadvection(dt, stepper, dev; kwargs...) + test_1layerqg_nonlinearadvection(dt, stepper, dev::Device=CPU(); + add_topography = false, + add_background_flow = false, + background_flow_vary_in_y = true) -Tests the advection term in the SingleLayerQG module by timestepping a +Tests the advection term in the `SingleLayerQG` module by timestepping a test problem with timestep dt and timestepper identified by the string stepper. -The test problem is derived by picking a solution ζf (with associated -streamfunction ψf) for which the advection term J(ψf, ζf) is non-zero. Next, a -forcing Ff is derived according to Ff = ∂ζf/∂t + J(ψf, ζf) - ν∇²ζf. One solution -to the vorticity equation forced by this Ff is then ζf. (This solution may not +The test problem is derived by picking a solution qf (with associated +streamfunction ψf) for which the advection term J(qf, ζf) is non-zero. Next, a +forcing Ff is derived according to Ff = ∂qf/∂t + J(ψf, qf) - ν∇²qf. One solution +to the vorticity equation forced by this Ff is then qf. (This solution may not be realized, at least at long times, if it is unstable.) + +We can optionally add an imposed mean flow `U(y)` via the kwarg `add_background_flow = true` +or add topography via kwarg `add_topography = true`. """ -function test_1layerqg_advection(dt, stepper, dev::Device=CPU(); n=128, L=2π, ν=1e-2, nν=1, μ=0.0) +function test_1layerqg_nonlinearadvection(dt, stepper, dev::Device=CPU(); + add_topography = false, + add_background_flow = false, + background_flow_vary_in_y = true) n, L = 128, 2π ν, nν = 1e-2, 1 μ = 0.0 @@ -275,13 +286,42 @@ function test_1layerqg_advection(dt, stepper, dev::Device=CPU(); n=128, L=2π, grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) + if add_topography == true + η₀ = 0.4 + elseif add_topography == false + η₀ = 0 + end + + η(x, y, η₀) = η₀ * cos(10x) * cos(10y) + η_array = @. η(x, y, η₀) + ψf = @. sin(2x) * cos(2y) + 2sin(x) * cos(3y) qf = @. -8sin(2x) * cos(2y) - 20sin(x) * cos(3y) - Ff = @. -( - ν*( 64sin(2x) * cos(2y) + 200sin(x) * cos(3y) ) - + 8 * ( cos(x) * cos(3y) * sin(2x) * sin(2y) - 3cos(2x) * cos(2y) * sin(x) * sin(3y) ) - ) + # F = J(ψ, q + η) - ν ∇²q + U ∂(q+η)/∂x - v ∂²U/∂y² + # where q = ∇²ψ + Ff = @. (- ν * ( 64sin(2x) * cos(2y) + 200sin(x) * cos(3y) ) + - 8 * (cos(x) * cos(3y) * sin(2x) * sin(2y) - 3cos(2x) * cos(2y) * sin(x) * sin(3y)) + - 20η₀ * (sin(10x) * cos(10y) * (sin(2x) * sin(2y) + 3sin(x) * sin(3y)) + + cos(10x) * sin(10y) * (cos(2x) * cos(2y) + cos(x) * cos(3y)))) + + U_amplitude = 0.2 + + U = 0 + Uyy = 0 + + if add_background_flow == true && background_flow_vary_in_y == false + U = U_amplitude + elseif add_background_flow == true && background_flow_vary_in_y == true + U = @. U_amplitude / cosh(2y)^2 + U = device_array(dev)(U) + Uh = rfft(U) + Uyy = irfft(- grid.l.^2 .* Uh, grid.nx) # ∂²U/∂y² + end + + @. Ff += (- U * (16cos(2x) * cos(2y) + 20cos(x) * cos(3y)) + - Uyy * (2cos(2x) * cos(2y) + 2cos(x) * cos(3y)) + - η₀ * U * 10sin(10x) * cos(10y)) Ffh = rfft(Ff) @@ -290,14 +330,71 @@ function test_1layerqg_advection(dt, stepper, dev::Device=CPU(); n=128, L=2π, return nothing end - prob = SingleLayerQG.Problem(dev; nx=n, Lx=L, ν=ν, nν=nν, μ=μ, dt=dt, stepper=stepper, calcF=calcF!) + if add_background_flow == false + U₀ = 0 + elseif add_background_flow == true && background_flow_vary_in_y == true + U₀ = U[1, :] + elseif add_background_flow == true && background_flow_vary_in_y == false + U₀ = U_amplitude + end + + prob = SingleLayerQG.Problem(dev; nx=n, Lx=L, U = U₀, eta = η_array, ν=ν, nν=nν, μ=μ, dt=dt, stepper=stepper, calcF=calcF!) SingleLayerQG.set_q!(prob, qf) stepforward!(prob, round(Int, nt)) SingleLayerQG.updatevars!(prob) - + + return isapprox(prob.vars.q, qf, rtol=rtol_singlelayerqg) +end + +""" + test_1layerqg_nonlinearadvection_deformation(dt, stepper, dev::Device=CPU(); + deformation_radius=1.23) + +Same as `test_1layerqg_nonlinearadvection` but with finite deformation radius. +""" +function test_1layerqg_nonlinearadvection_deformation(dt, stepper, dev::Device=CPU(); + deformation_radius=1.23) + n, L = 128, 2π + ν, nν = 1e-2, 1 + μ = 0.0 + tf = 1.0 + nt = round(Int, tf/dt) + + grid = TwoDGrid(dev; nx=n, Lx=L) + k₀, l₀ = 2π / grid.Lx, 2π / grid.Ly # fundamental wavenumbers + x, y = gridpoints(grid) + + η₀ = 0.4 + η(x, y) = η₀ * cos(10x) * cos(10y) + ψf = @. sin(2x) * cos(2y) + 2sin(x) * cos(3y) + qf = @. - (2^2 + 2^2) * sin(2x) * cos(2y) - (1^2 + 3^2) * 2sin(x) * cos(3y) - 1/deformation_radius^2 * ψf + + # F = J(ψ, q + η) - ν ∇²q + # where q = ∇²ψ - ψ/ℓ² + Ff = @. (- ν * (64sin(2x) * cos(2y) + 200sin(x) * cos(3y) + + (8sin(2x) * cos(2y) + 20sin(x) * cos(3y)) / deformation_radius^2) + - 8 * (cos(x) * cos(3y) * sin(2x) * sin(2y) - 3cos(2x) * cos(2y) * sin(x) * sin(3y)) + - 20η₀ * (cos(10y) * sin(10x) * (sin(2x) * sin(2y) + 3sin(x) * sin(3y)) + + cos(10x) * sin(10y) * (cos(2x) * cos(2y) + cos(x) * cos(3y)))) + + Ffh = rfft(Ff) + + function calcF!(Fh, sol, t, clock, vars, params, grid) + Fh .= Ffh + return nothing + end + + prob = SingleLayerQG.Problem(dev; nx=n, Lx=L, eta=η, ν=ν, nν=nν, μ=μ, dt=dt, deformation_radius=deformation_radius, stepper=stepper, calcF=calcF!) + + SingleLayerQG.set_q!(prob, qf) + + stepforward!(prob, round(Int, nt)) + + SingleLayerQG.updatevars!(prob) + return isapprox(prob.vars.q, qf, rtol=rtol_singlelayerqg) end @@ -309,50 +406,59 @@ Tests the energy and enstrophy function for a SingleLayerQG problem. function test_1layerqg_energyenstrophy_BarotropicQG(dev::Device=CPU()) nx, Lx = 64, 2π ny, Ly = 64, 3π - grid = TwoDGrid(dev; nx, Lx, ny, Ly) - k₀, l₀ = 2π/Lx, 2π/Ly # fundamental wavenumbers + + prob = SingleLayerQG.Problem(dev; nx=nx, Lx=Lx, ny=ny, Ly=Ly, stepper="ForwardEuler") + grid = prob.grid + + k₀, l₀ = 2π / grid.Lx, 2π / grid.Ly # fundamental wavenumbers x, y = gridpoints(grid) energy_calc = 29/9 enstrophy_calc = 10885/648 - η = @. cos(10k₀ * x) * cos(10l₀ * y) + η(x, y) = cos(10k₀ * x) * cos(10l₀ * y) ψ₀ = @. sin(2k₀ * x) * cos(2l₀ * y) + 2sin(k₀ * x) * cos(3l₀ * y) q₀ = @. - ((2k₀)^2 + (2l₀)^2) * sin(2k₀ * x) * cos(2l₀ * y) - (k₀^2 + (3l₀)^2) * 2sin(k₀ * x) * cos(3l₀*y) prob = SingleLayerQG.Problem(dev; nx=nx, Lx=Lx, ny=ny, Ly=Ly, eta=η, stepper="ForwardEuler") sol, clock, vars, params, grid = prob.sol, prob.clock, prob.vars, prob.params, prob.grid + SingleLayerQG.set_q!(prob, q₀) SingleLayerQG.updatevars!(prob) energyq₀ = SingleLayerQG.energy(prob) enstrophyq₀ = SingleLayerQG.enstrophy(prob) - return isapprox(energyq₀, energy_calc, rtol=rtol_singlelayerqg) && isapprox(enstrophyq₀, enstrophy_calc, rtol=rtol_singlelayerqg) && SingleLayerQG.potential_energy(prob)==0 && - SingleLayerQG.addforcing!(prob.timestepper.N, sol, clock.t, clock, vars, params, grid) == nothing + return (isapprox(energyq₀, energy_calc, rtol=rtol_singlelayerqg) && + isapprox(enstrophyq₀, enstrophy_calc, rtol=rtol_singlelayerqg) && + SingleLayerQG.potential_energy(prob)==0 && + SingleLayerQG.addforcing!(prob.timestepper.N, sol, clock.t, clock, vars, params, grid) === nothing) end """ - test_1layerqg_energies_EquivalentBarotropicQG(dev) + test_1layerqg_energies_EquivalentBarotropicQG(dev; deformation_radius=1.23) Tests the kinetic and potential energy for an equivalent barotropic SingleLayerQG problem. """ function test_1layerqg_energies_EquivalentBarotropicQG(dev; deformation_radius=1.23) nx, Lx = 64, 2π ny, Ly = 64, 3π - grid = TwoDGrid(dev; nx, Lx, ny, Ly) - k₀, l₀ = 2π/Lx, 2π/Ly # fundamental wavenumbers + + prob = SingleLayerQG.Problem(dev; nx, Lx, ny, Ly, deformation_radius, stepper="ForwardEuler") + grid = prob.grid + + k₀, l₀ = 2π / grid.Lx, 2π / grid.Ly # fundamental wavenumbers x, y = gridpoints(grid) kinetic_energy_calc = 29/9 potential_energy_calc = 5/(8*deformation_radius^2) energy_calc = kinetic_energy_calc + potential_energy_calc - - η = @. cos(10k₀ * x) * cos(10l₀ * y) + + η(x, y) = cos(10k₀ * x) * cos(10l₀ * y) ψ₀ = @. sin(2k₀ * x) * cos(2l₀ * y) + 2sin(k₀ * x) * cos(3l₀ * y) q₀ = @. - ((2k₀)^2 + (2l₀)^2) * sin(2k₀ * x) * cos(2l₀ * y) - (k₀^2 + (3l₀)^2) * 2sin(k₀ * x) * cos(3l₀*y) - 1/deformation_radius^2 * ψ₀ - prob = SingleLayerQG.Problem(dev; nx=nx, Lx=Lx, ny=ny, Ly=Ly, eta=η, deformation_radius=deformation_radius, stepper="ForwardEuler") + prob = SingleLayerQG.Problem(dev; nx, Lx, ny, Ly, eta=η, deformation_radius, stepper="ForwardEuler") sol, clock, vars, params, grid = prob.sol, prob.clock, prob.vars, prob.params, prob.grid SingleLayerQG.set_q!(prob, q₀) SingleLayerQG.updatevars!(prob) @@ -361,21 +467,24 @@ function test_1layerqg_energies_EquivalentBarotropicQG(dev; deformation_radius=1 potential_energyq₀ = SingleLayerQG.potential_energy(prob) energyq₀ = SingleLayerQG.energy(prob) - return isapprox(kinetic_energyq₀, kinetic_energy_calc, rtol=rtol_singlelayerqg) && isapprox(potential_energyq₀, potential_energy_calc, rtol=rtol_singlelayerqg) && isapprox(energyq₀, energy_calc, rtol=rtol_singlelayerqg) && - SingleLayerQG.addforcing!(prob.timestepper.N, sol, clock.t, clock, vars, params, grid) == nothing + return (isapprox(kinetic_energyq₀, kinetic_energy_calc, rtol=rtol_singlelayerqg) && + isapprox(potential_energyq₀, potential_energy_calc, rtol=rtol_singlelayerqg) && + isapprox(energyq₀, energy_calc, rtol=rtol_singlelayerqg) && + SingleLayerQG.addforcing!(prob.timestepper.N, sol, clock.t, clock, vars, params, grid) === nothing) end """ - test_1layerqg_problemtype(dev, T; deformation_radius=Inf) + test_1layerqg_problemtype(dev, T; deformation_radius=Inf, U=0) -Tests the SingleLayerQG problem constructor for different DataType `T`. +Test the SingleLayerQG problem constructor for different DataType `T`. """ -function test_1layerqg_problemtype(dev, T; deformation_radius=Inf) - prob = SingleLayerQG.Problem(dev; T=T, deformation_radius=deformation_radius) +function test_1layerqg_problemtype(dev, T; deformation_radius=Inf, U=0) + prob = SingleLayerQG.Problem(dev; T, deformation_radius, U) A = device_array(dev) - - return (typeof(prob.sol)<:A{Complex{T}, 2} && typeof(prob.grid.Lx)==T && eltype(prob.grid.x)==T && typeof(prob.vars.u)<:A{T, 2}) + + return (typeof(prob.sol)<:A{Complex{T}, 2} && typeof(prob.grid.Lx)==T && + eltype(prob.grid.x)==T && typeof(prob.vars.u)<:A{T, 2} && typeof(prob.params.U)==T) end function test_streamfunctionfrompv(dev; deformation_radius=1.23) @@ -395,7 +504,6 @@ function test_streamfunctionfrompv(dev; deformation_radius=1.23) SingleLayerQG.set_q!(prob_equivalentbarotropicQG, q_equivalentbarotropic) SingleLayerQG.streamfunctionfrompv!(prob_barotropicQG.vars.ψh, prob_barotropicQG.vars.qh, prob_barotropicQG.params, prob_barotropicQG.grid) - SingleLayerQG.streamfunctionfrompv!(prob_equivalentbarotropicQG.vars.ψh, prob_equivalentbarotropicQG.vars.qh, prob_equivalentbarotropicQG.params, prob_equivalentbarotropicQG.grid) return (prob_barotropicQG.vars.ψ ≈ ψ && prob_equivalentbarotropicQG.vars.ψ ≈ ψ) @@ -435,4 +543,4 @@ function test_1layerqg_enstrophy_drag(dev; deformation_radius=2.23) prob = SingleLayerQG.Problem(dev; nx=16, deformation_radius=deformation_radius) SingleLayerQG.enstrophy_drag(prob) return nothing -end \ No newline at end of file +end From 1dd98e4b9d7e6fe421d1423db8bad70fe646d07f Mon Sep 17 00:00:00 2001 From: "Gregory L. Wagner" Date: Sat, 7 Sep 2024 15:09:28 -0600 Subject: [PATCH 3/7] Clarify `examples/README` (#366) I thought it was confusing that the directory is `examples/` but the README says to navigate to `examples/directory` --- examples/README.md | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/examples/README.md b/examples/README.md index de17ba17..71371a5b 100644 --- a/examples/README.md +++ b/examples/README.md @@ -1,9 +1,7 @@ # GeophysicalFlows.jl/examples - These are some basic examples of the various modules included in GeophysicalFlows.jl. The best way to go through the examples by browsing them within the package's documentation . - ## Run the examples You can run the examples in an executable environment via [Binder](https://mybinder.org) by clicking on the [![badge](https://img.shields.io/badge/binder-badge-579ACA.svg?logo=data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAFkAAABZCAMAAABi1XidAAAB8lBMVEX///9XmsrmZYH1olJXmsr1olJXmsrmZYH1olJXmsr1olJXmsrmZYH1olL1olJXmsr1olJXmsrmZYH1olL1olJXmsrmZYH1olJXmsr1olL1olJXmsrmZYH1olL1olJXmsrmZYH1olL1olL0nFf1olJXmsrmZYH1olJXmsq8dZb1olJXmsrmZYH1olJXmspXmspXmsr1olL1olJXmsrmZYH1olJXmsr1olL1olJXmsrmZYH1olL1olLeaIVXmsrmZYH1olL1olL1olJXmsrmZYH1olLna31Xmsr1olJXmsr1olJXmsrmZYH1olLqoVr1olJXmsr1olJXmsrmZYH1olL1olKkfaPobXvviGabgadXmsqThKuofKHmZ4Dobnr1olJXmsr1olJXmspXmsr1olJXmsrfZ4TuhWn1olL1olJXmsqBi7X1olJXmspZmslbmMhbmsdemsVfl8ZgmsNim8Jpk8F0m7R4m7F5nLB6jbh7jbiDirOEibOGnKaMhq+PnaCVg6qWg6qegKaff6WhnpKofKGtnomxeZy3noG6dZi+n3vCcpPDcpPGn3bLb4/Mb47UbIrVa4rYoGjdaIbeaIXhoWHmZYHobXvpcHjqdHXreHLroVrsfG/uhGnuh2bwj2Hxk17yl1vzmljzm1j0nlX1olL3AJXWAAAAbXRSTlMAEBAQHx8gICAuLjAwMDw9PUBAQEpQUFBXV1hgYGBkcHBwcXl8gICAgoiIkJCQlJicnJ2goKCmqK+wsLC4usDAwMjP0NDQ1NbW3Nzg4ODi5+3v8PDw8/T09PX29vb39/f5+fr7+/z8/Pz9/v7+zczCxgAABC5JREFUeAHN1ul3k0UUBvCb1CTVpmpaitAGSLSpSuKCLWpbTKNJFGlcSMAFF63iUmRccNG6gLbuxkXU66JAUef/9LSpmXnyLr3T5AO/rzl5zj137p136BISy44fKJXuGN/d19PUfYeO67Znqtf2KH33Id1psXoFdW30sPZ1sMvs2D060AHqws4FHeJojLZqnw53cmfvg+XR8mC0OEjuxrXEkX5ydeVJLVIlV0e10PXk5k7dYeHu7Cj1j+49uKg7uLU61tGLw1lq27ugQYlclHC4bgv7VQ+TAyj5Zc/UjsPvs1sd5cWryWObtvWT2EPa4rtnWW3JkpjggEpbOsPr7F7EyNewtpBIslA7p43HCsnwooXTEc3UmPmCNn5lrqTJxy6nRmcavGZVt/3Da2pD5NHvsOHJCrdc1G2r3DITpU7yic7w/7Rxnjc0kt5GC4djiv2Sz3Fb2iEZg41/ddsFDoyuYrIkmFehz0HR2thPgQqMyQYb2OtB0WxsZ3BeG3+wpRb1vzl2UYBog8FfGhttFKjtAclnZYrRo9ryG9uG/FZQU4AEg8ZE9LjGMzTmqKXPLnlWVnIlQQTvxJf8ip7VgjZjyVPrjw1te5otM7RmP7xm+sK2Gv9I8Gi++BRbEkR9EBw8zRUcKxwp73xkaLiqQb+kGduJTNHG72zcW9LoJgqQxpP3/Tj//c3yB0tqzaml05/+orHLksVO+95kX7/7qgJvnjlrfr2Ggsyx0eoy9uPzN5SPd86aXggOsEKW2Prz7du3VID3/tzs/sSRs2w7ovVHKtjrX2pd7ZMlTxAYfBAL9jiDwfLkq55Tm7ifhMlTGPyCAs7RFRhn47JnlcB9RM5T97ASuZXIcVNuUDIndpDbdsfrqsOppeXl5Y+XVKdjFCTh+zGaVuj0d9zy05PPK3QzBamxdwtTCrzyg/2Rvf2EstUjordGwa/kx9mSJLr8mLLtCW8HHGJc2R5hS219IiF6PnTusOqcMl57gm0Z8kanKMAQg0qSyuZfn7zItsbGyO9QlnxY0eCuD1XL2ys/MsrQhltE7Ug0uFOzufJFE2PxBo/YAx8XPPdDwWN0MrDRYIZF0mSMKCNHgaIVFoBbNoLJ7tEQDKxGF0kcLQimojCZopv0OkNOyWCCg9XMVAi7ARJzQdM2QUh0gmBozjc3Skg6dSBRqDGYSUOu66Zg+I2fNZs/M3/f/Grl/XnyF1Gw3VKCez0PN5IUfFLqvgUN4C0qNqYs5YhPL+aVZYDE4IpUk57oSFnJm4FyCqqOE0jhY2SMyLFoo56zyo6becOS5UVDdj7Vih0zp+tcMhwRpBeLyqtIjlJKAIZSbI8SGSF3k0pA3mR5tHuwPFoa7N7reoq2bqCsAk1HqCu5uvI1n6JuRXI+S1Mco54YmYTwcn6Aeic+kssXi8XpXC4V3t7/ADuTNKaQJdScAAAAAElFTkSuQmCC)](https://mybinder.org) @@ -12,7 +10,7 @@ You can run the examples in an executable environment via [Binder](https://mybin Alternatively, you can run these scripts directly using the `.toml` files in this directory. To do that, first open julia and activate this directory's project, e.g., ``` -$ julia --project="path/to/examples/directory" +$ julia --project="path/to/examples" ``` Then instantiate the project in this directory, i.e., ```julia @@ -20,5 +18,5 @@ Then instantiate the project in this directory, i.e., ``` to install dependencies. Then run any of the examples by ```julia -include("path/to/examples/directory/example_script.jl") +include("path/to/examples/example_script.jl") ``` From f0714b1cf241453d1ea9dbedc73246b06bf66373 Mon Sep 17 00:00:00 2001 From: "Navid C. Constantinou" Date: Tue, 22 Oct 2024 14:30:11 +1100 Subject: [PATCH 4/7] Julia v1.11 for CI + update installation Julia version recommendations (#371) * julia v1.11 for ci * cancel-workflow-action@0.12 * cancel-workflow-action@0.12.1 * remove random space * add windows CI via gh actions * remove appveyor * add docs * remove Documenter action * dont test nightly on win * update compathelper * add julia v1.10 for macOS * add julia v1.10 for macOS * tweak ci * tweak ci * update installation recommendation for julia versions * tweak ci --- .buildkite/pipeline.yml | 12 +++--- .github/workflows/CI.yml | 54 +++++++++++++++++++++++++++ .github/workflows/CompatHelper.yml | 44 +++++++++++++++++----- .github/workflows/Documenter.yml | 29 -------------- README.md | 10 ++--- appveyor.yml | 36 ------------------ docs/src/installation_instructions.md | 12 +++--- 7 files changed, 105 insertions(+), 92 deletions(-) create mode 100644 .github/workflows/CI.yml delete mode 100644 .github/workflows/Documenter.yml delete mode 100644 appveyor.yml diff --git a/.buildkite/pipeline.yml b/.buildkite/pipeline.yml index e938b001..e7f4d573 100644 --- a/.buildkite/pipeline.yml +++ b/.buildkite/pipeline.yml @@ -3,24 +3,24 @@ env: SECRET_CODECOV_TOKEN: "Ak2mVTxXnhkPNc096ImDdp7bOc4zGNTqFEEDaGMwAgYPr28g5dyMbslh8B/ad4NQHXVL1MXQ3zrUfGgMBRq+mmqRaAe13FI4Go9uCas6bzdZXE23ExiLzBmqVRNRf8GqEcpGL7BBreohC0cnfI0SVMiIJDCJXX9YsXJtlcpk1glQFMEFI5V6cpFe9K2l5xoUNQ4179ZYoJUMAy/aylQx/UdQuw527FjHQUsi5/dFtWzMqeys0secNa9alLvJCQdIX9OqPjmBYvuIIVXCR7vlZoH8PgXwEj7wbdp8/V31+wlLQI9WePcsJxoOybtLTOlwwfw4jWLAttDZYqnqiLVp3Q==;U2FsdGVkX18sNLCManU1B/jI5kh4LhSi69MFXljHSp9yWrN7u5d196K/XrELwb8ksbamyKeHjIvDIopwD55dbw==" steps: - - label: "🥝 Julia 1.6" + - label: "🥑 Julia 1.10" plugins: - JuliaCI/julia#v1: - version: '1.6' + version: '1.10' - JuliaCI/julia-test#v1: ~ - - JuliaCI/julia-coverage#v1: - codecov: true agents: queue: "juliagpu" cuda: "*" if: build.message !~ /\[skip tests\]/ timeout_in_minutes: 60 - - label: "🥑 Julia 1.10" + - label: "🥝 Julia 1.11" plugins: - JuliaCI/julia#v1: - version: '1.10' + version: '1.11' - JuliaCI/julia-test#v1: ~ + - JuliaCI/julia-coverage#v1: + codecov: true agents: queue: "juliagpu" cuda: "*" diff --git a/.github/workflows/CI.yml b/.github/workflows/CI.yml new file mode 100644 index 00000000..587cf11a --- /dev/null +++ b/.github/workflows/CI.yml @@ -0,0 +1,54 @@ +name: CI +on: + push: + branches: [main] + tags: ["*"] + pull_request: +jobs: + test: + name: Julia ${{ matrix.version }} - ${{ matrix.os }} - ${{ matrix.arch }} - ${{ github.event_name }} + runs-on: ${{ matrix.os }} + strategy: + fail-fast: false + matrix: + version: + - "1.10" # LTS + - "1.11" + os: + - macOS-latest + - windows-latest + arch: + - x64 + - aarch64 + - arm64 + exclude: + # Test 32-bit only on Linux + - os: macOS-latest + arch: x64 + - os: windows-latest + arch: aarch64 + - os: windows-latest + arch: arm64 + steps: + - uses: actions/checkout@v4 + - uses: julia-actions/setup-julia@v2 + with: + version: ${{ matrix.version }} + arch: ${{ matrix.arch }} + - uses: julia-actions/cache@v2 + - uses: julia-actions/julia-buildpkg@v1 + - uses: julia-actions/julia-runtest@v1 + docs: + name: Documentation + runs-on: ubuntu-latest + steps: + - uses: actions/checkout@v4 + - uses: julia-actions/setup-julia@v2 + with: + version: '1.11' + - uses: julia-actions/julia-buildpkg@latest + - uses: julia-actions/julia-docdeploy@latest + env: + GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }} + DOCUMENTER_KEY: ${{ secrets.DOCUMENTER_KEY }} + JULIA_DEBUG: Documenter diff --git a/.github/workflows/CompatHelper.yml b/.github/workflows/CompatHelper.yml index 821653f2..c12d0638 100644 --- a/.github/workflows/CompatHelper.yml +++ b/.github/workflows/CompatHelper.yml @@ -1,19 +1,45 @@ name: CompatHelper - on: schedule: - - cron: '00 00 * * *' - + - cron: 0 0 * * * + workflow_dispatch: +permissions: + contents: write + pull-requests: write jobs: CompatHelper: runs-on: ubuntu-latest steps: - - uses: julia-actions/setup-julia@latest + - name: Check if Julia is already available in the PATH + id: julia_in_path + run: which julia + continue-on-error: true + - name: Install Julia, but only if it is not already available in the PATH + uses: julia-actions/setup-julia@v2 with: - version: '1.10' - - name: Pkg.add("CompatHelper") - run: julia -e 'using Pkg; Pkg.add("CompatHelper")' - - name: CompatHelper.main() + version: '1' + # arch: ${{ runner.arch }} + if: steps.julia_in_path.outcome != 'success' + - name: "Add the General registry via Git" + run: | + import Pkg + ENV["JULIA_PKG_SERVER"] = "" + Pkg.Registry.add("General") + shell: julia --color=yes {0} + - name: "Install CompatHelper" + run: | + import Pkg + name = "CompatHelper" + uuid = "aa819f21-2bde-4658-8897-bab36330d9b7" + version = "3" + Pkg.add(; name, uuid, version) + shell: julia --color=yes {0} + - name: "Run CompatHelper" + run: | + import CompatHelper + CompatHelper.main() + shell: julia --color=yes {0} env: GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }} - run: julia -e 'using CompatHelper; CompatHelper.main()' + COMPATHELPER_PRIV: ${{ secrets.DOCUMENTER_KEY }} + # COMPATHELPER_PRIV: ${{ secrets.COMPATHELPER_PRIV }} diff --git a/.github/workflows/Documenter.yml b/.github/workflows/Documenter.yml deleted file mode 100644 index d85b18b9..00000000 --- a/.github/workflows/Documenter.yml +++ /dev/null @@ -1,29 +0,0 @@ -name: Documentation - -on: - push: - branches: - - main - tags: '*' - pull_request: - -jobs: - build: - runs-on: ubuntu-latest - steps: - - name: Cancel Previous Runs - uses: styfle/cancel-workflow-action@0.11.0 - with: - access_token: ${{ github.token }} - - uses: actions/checkout@v4 - - uses: julia-actions/setup-julia@latest - with: - version: '1.10' - - name: Install dependencies - run: julia --project=docs/ -e 'using Pkg; Pkg.develop(PackageSpec(path=pwd())); Pkg.instantiate()' - - name: Build and deploy - env: - GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }} # For authentication with GitHub Actions token - DOCUMENTER_KEY: ${{ secrets.DOCUMENTER_KEY }} # For authentication with SSH deploy key - JULIA_DEBUG: Documenter - run: julia --color=yes --project=docs/ docs/make.jl diff --git a/README.md b/README.md index 48b40e04..ec906796 100644 --- a/README.md +++ b/README.md @@ -43,14 +43,12 @@ To install, use Julia's built-in package manager (accessed by pressing `]` in t ```julia julia> ] -(v1.6) pkg> add GeophysicalFlows -(v1.6) pkg> instantiate +(v1.10) pkg> add GeophysicalFlows +(v1.10) pkg> instantiate ``` -The most recent version of GeophysicalFlows.jl requires Julia v1.6 (the current long-term-release) or later. _We strongly urge you to use this version._ - -The latest version that is compatible with Julia v1.5 is GeophysicalFlows.jl v0.13.1. - +GeophysicalFlows.jl requires Julia v1.6 or later. However, the package has continuous integration testing on +Julia v1.10 (the current long-term release) and v1.11. _We strongly urge you to use one of these Julia versions._ ## Examples diff --git a/appveyor.yml b/appveyor.yml deleted file mode 100644 index 1dc6d083..00000000 --- a/appveyor.yml +++ /dev/null @@ -1,36 +0,0 @@ -environment: - matrix: - - julia_version: 1.6 - - julia_version: 1.10 - - julia_version: nightly - -platform: - - x64 # 64-bit - -# Uncomment the following lines to allow failures on nightly julia -# (tests will run but not make your overall status red) -matrix: - allow_failures: - - julia_version: nightly - -branches: - only: - - main - - /release-.*/ - -notifications: - - provider: Email - on_build_success: false - on_build_failure: false - on_build_status_changed: false - -install: - - ps: iex ((new-object net.webclient).DownloadString("https://raw.githubusercontent.com/JuliaCI/Appveyor.jl/version-1/bin/install.ps1")) - -build_script: - - echo "%JL_BUILD_SCRIPT%" - - C:\julia\bin\julia -e "%JL_BUILD_SCRIPT%" - -test_script: - - echo "%JL_TEST_SCRIPT%" - - C:\julia\bin\julia -e "%JL_TEST_SCRIPT%" diff --git a/docs/src/installation_instructions.md b/docs/src/installation_instructions.md index 48a6ed85..1fec893f 100644 --- a/docs/src/installation_instructions.md +++ b/docs/src/installation_instructions.md @@ -6,8 +6,8 @@ instantiate/build all the required dependencies ```julia julia> ] -(v1.6) pkg> add GeophysicalFlows -(v1.6) pkg> instantiate +(v1.10) pkg> add GeophysicalFlows +(v1.10) pkg> instantiate ``` We recommend installing GeophysicalFlows.jl with the built-in Julia package manager, because @@ -15,7 +15,7 @@ this installs a stable, tagged release. Later on, you can update GeophysicalFlow latest tagged release again via the package manager by typing ```julia -(v1.6) pkg> update GeophysicalFlows +(v1.10) pkg> update GeophysicalFlows ``` Note that some releases might induce breaking changes to certain modules. If after anything @@ -23,10 +23,10 @@ happens or your code stops working, please open an [issue](https://github.com/Fo or start a [discussion](https://github.com/FourierFlows/GeophysicalFlows.jl/discussions). We're more than happy to help with getting your simulations up and running. -!!! warn "Use Julia 1.6 or newer" +!!! warn "Julia 1.6 or newer required; Julia 1.10 or newer strongly encouraged" The latest GeophysicalFlows.jl requires at least Julia v1.6 to run. Installing GeophysicalFlows with an older version of Julia will install an older version of GeophysicalFlows.jl (the latest version compatible with your version of Julia). - - Last version compatible with Julia v1.5: GeophysicalFlows.jl v0.13.1 + GeophysicalFlows.jl is continuously tested on Julia v1.10 (the current long-term release) and v1.11. + _We strongly urge you to use one of these Julia versions._ From 700c675ef53975d5f8cb5d62159848b55d38f5aa Mon Sep 17 00:00:00 2001 From: "Navid C. Constantinou" Date: Tue, 22 Oct 2024 14:51:15 +1100 Subject: [PATCH 5/7] update paper reference (#369) --- docs/src/index.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/src/index.md b/docs/src/index.md index 6abfdb27..a1a03e90 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -92,7 +92,7 @@ The bibtex entry for the paper is: 1. He, J. and Wang, Y. (2024) Multiple states of two-dimensional turbulence above topography. arXiv preprint arXiv:2405.10826, doi:[10.48550/arXiv.2405.10826](https://doi.org/10.48550/arXiv.2405.10826). -1. Parfenyev, V., Blumenau, M., and Nikitin, I. (2024) Enhancing capabilities of particle image/tracking velocimetry with physics-informed neural networks. arXiv preprint arXiv:2404.01193, doi:[10.48550/arXiv.2404.01193](https://doi.org/10.48550/arXiv.2404.01193). +1. Parfenyev, V., Blumenau, M., and Nikitin, I. (2024) Inferring parameters and reconstruction of two-dimensional turbulent flows with physics-informed neural networks. _Jetp Lett._, doi:[10.1134/S0021364024602203](https://doi.org/10.1134/S0021364024602203). 1. Shokar, I. J. S., Kerswell, R. R., and Haynes, P. H. (2024) Stochastic latent transformer: Efficient modeling of stochastically forced zonal jets. _Journal of Advances in Modeling Earth Systems_, **16**, e2023MS004177, doi:[10.1029/2023MS004177](https://doi.org/10.1029/2023MS004177). From 960518b0de417a44a449544f208a75d16c845f46 Mon Sep 17 00:00:00 2001 From: "Navid C. Constantinou" Date: Tue, 22 Oct 2024 14:51:33 +1100 Subject: [PATCH 6/7] add Crowe and Sutyrin (#370) --- docs/src/index.md | 2 ++ 1 file changed, 2 insertions(+) diff --git a/docs/src/index.md b/docs/src/index.md index a1a03e90..c8c10357 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -84,6 +84,8 @@ The bibtex entry for the paper is: ## Papers using `GeophysicalFlows.jl` +1. Crowe, M. N. and Sutyrin, G. G. (2024) Symmetry breaking of two-layer eastward propagating dipoles. arXiv preprint arXiv.2410.14402, doi:[10.48550/arXiv.2410.14402](https://doi.org/10.48550/arXiv.2410.14402). + 1. Lobo, M., Griffies, S. M., and Zhang, W. (2024) Vertical structure of baroclinic instability in a three-layer quasi-geostrophic model over a sloping bottom. _ESS Open Archive_, doi:[10.22541/essoar.172166435.57577370/v1](https://doi.org/10.22541/essoar.172166435.57577370/v1). 1. Pudig, M. and Smith, K. S. (2024) Baroclinic turbulence above rough topography: The vortex gas and topographic turbulence regimes. _ESS Open Archive_, doi:[10.22541/essoar.171995116.60993353/v1](https://doi.org/10.22541/essoar.171995116.60993353/v1). From 8f9438c6c3e572462dd5cff96b2acb905f03dae5 Mon Sep 17 00:00:00 2001 From: "Navid C. Constantinou" Date: Tue, 22 Oct 2024 17:01:10 +1100 Subject: [PATCH 7/7] Drop appveyor badge (#372) * Drop appveyor badge * add CI badge alt label --- README.md | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/README.md b/README.md index ec906796..e27b777b 100644 --- a/README.md +++ b/README.md @@ -10,8 +10,8 @@ Buildkite CPU+GPU build status - - Build Status for Window + + CI Status stable docs