From 0676329c6aed72dda5f4ab441b717d4461aaef09 Mon Sep 17 00:00:00 2001 From: "Navid C. Constantinou" Date: Fri, 2 Sep 2022 13:46:29 +1000 Subject: [PATCH 1/7] update to latest grid refactor --- Project.toml | 2 +- docs/src/aliasing.md | 2 +- examples/barotropicqgql_betaforced.jl | 4 +- examples/multilayerqg_2layer.jl | 4 +- examples/singlelayerqg_betadecay.jl | 4 +- examples/singlelayerqg_betaforced.jl | 4 +- examples/singlelayerqg_decaying_topography.jl | 4 +- .../twodnavierstokes_stochasticforcing.jl | 4 +- ...dnavierstokes_stochasticforcing_budgets.jl | 4 +- src/barotropicqgql.jl | 42 +++++++------- src/multilayerqg.jl | 56 ++++++++++--------- test/test_barotropicqgql.jl | 6 +- test/test_multilayerqg.jl | 6 +- test/test_singlelayerqg.jl | 6 +- test/test_surfaceqg.jl | 12 ++-- test/test_twodnavierstokes.jl | 14 ++--- 16 files changed, 92 insertions(+), 82 deletions(-) diff --git a/Project.toml b/Project.toml index f0b6ff2e..31715bbd 100644 --- a/Project.toml +++ b/Project.toml @@ -21,7 +21,7 @@ Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" CUDA = "^1, ^2.4.2, 3.0.0 - 3.6.4, ^3.7.1" DocStringExtensions = "^0.8, 0.9" FFTW = "^1" -FourierFlows = "^0.9.4" +FourierFlows = "^0.10.0" JLD2 = "^0.1, ^0.2, ^0.3, ^0.4" Reexport = "^0.2, ^1" SpecialFunctions = "^0.10, ^1, 2" diff --git a/docs/src/aliasing.md b/docs/src/aliasing.md index fac8a21b..11162558 100644 --- a/docs/src/aliasing.md +++ b/docs/src/aliasing.md @@ -11,7 +11,7 @@ in Fourier space before transforming to physical space to compute nonlinear term Users can construct a `grid` with different `aliased_fraction` via ```julia - julia> grid = OneDGrid(64, 2π, aliased_fraction=1/2) + julia> grid = OneDGrid(; nx=64, Lx=2π, aliased_fraction=1/2) julia> OneDimensionalGrid ├─────────── Device: CPU diff --git a/examples/barotropicqgql_betaforced.jl b/examples/barotropicqgql_betaforced.jl index 23765e51..0b053f35 100644 --- a/examples/barotropicqgql_betaforced.jl +++ b/examples/barotropicqgql_betaforced.jl @@ -62,7 +62,7 @@ forcing_wavenumber = 14.0 * 2π/L # the forcing wavenumber, `k_f`, for a spectr forcing_bandwidth = 1.5 * 2π/L # the width of the forcing spectrum, `δ_f` ε = 0.001 # energy input rate by the forcing -grid = TwoDGrid(dev, n, L) +grid = TwoDGrid(dev; nx=n, Lx=L) K = @. sqrt(grid.Krsq) # a 2D array with the total wavenumber @@ -133,7 +133,7 @@ fig # ## Setting initial conditions # Our initial condition is simply fluid at rest. -BarotropicQGQL.set_zeta!(prob, ArrayType(dev)(zeros(grid.nx, grid.ny))) +BarotropicQGQL.set_zeta!(prob, device_array(dev)(zeros(grid.nx, grid.ny))) nothing # hide # ## Diagnostics diff --git a/examples/multilayerqg_2layer.jl b/examples/multilayerqg_2layer.jl index 68caa332..e79edfa4 100644 --- a/examples/multilayerqg_2layer.jl +++ b/examples/multilayerqg_2layer.jl @@ -76,11 +76,11 @@ nothing # hide # Our initial condition is some small-amplitude random noise. We smooth our initial # condidtion using the `timestepper`'s high-wavenumber `filter`. # -# `ArrayType()` function returns the array type appropriate for the device, i.e., `Array` for +# `device_array()` function returns the array type appropriate for the device, i.e., `Array` for # `dev = CPU()` and `CuArray` for `dev = GPU()`. seed!(1234) # reset of the random number generator for reproducibility -q₀ = 1e-2 * ArrayType(dev)(randn((grid.nx, grid.ny, nlayers))) +q₀ = 1e-2 * device_array(dev)(randn((grid.nx, grid.ny, nlayers))) q₀h = prob.timestepper.filter .* rfft(q₀, (1, 2)) # apply rfft only in dims=1, 2 q₀ = irfft(q₀h, grid.nx, (1, 2)) # apply irfft only in dims=1, 2 diff --git a/examples/singlelayerqg_betadecay.jl b/examples/singlelayerqg_betadecay.jl index 3d275ce0..9e8a1c48 100644 --- a/examples/singlelayerqg_betadecay.jl +++ b/examples/singlelayerqg_betadecay.jl @@ -64,7 +64,7 @@ nothing # hide # Our initial condition consist of a flow that has power only at wavenumbers with # ``6 < \frac{L}{2\pi} \sqrt{k_x^2 + k_y^2} < 10`` and initial energy ``E_0``. -# `ArrayType()` function returns the array type appropriate for the device, i.e., `Array` for +# `device_array()` function returns the array type appropriate for the device, i.e., `Array` for # `dev = CPU()` and `CuArray` for `dev = GPU()`. E₀ = 0.08 # energy of initial condition @@ -72,7 +72,7 @@ E₀ = 0.08 # energy of initial condition K = @. sqrt(grid.Krsq) # a 2D array with the total wavenumber Random.seed!(1234) -q₀h = ArrayType(dev)(randn(Complex{eltype(grid)}, size(sol))) +q₀h = device_array(dev)(randn(Complex{eltype(grid)}, size(sol))) @. q₀h = ifelse(K < 6 * 2π/L, 0, q₀h) @. q₀h = ifelse(K > 10 * 2π/L, 0, q₀h) @. q₀h[1, :] = 0 # remove any power from zonal wavenumber k=0 diff --git a/examples/singlelayerqg_betaforced.jl b/examples/singlelayerqg_betaforced.jl index de8ea65d..4337df5c 100644 --- a/examples/singlelayerqg_betaforced.jl +++ b/examples/singlelayerqg_betaforced.jl @@ -63,7 +63,7 @@ forcing_wavenumber = 14.0 * 2π/L # the forcing wavenumber, `k_f`, for a spectr forcing_bandwidth = 1.5 * 2π/L # the width of the forcing spectrum, `δ_f` ε = 0.001 # energy input rate by the forcing -grid = TwoDGrid(dev, n, L) +grid = TwoDGrid(dev; nx=n, Lx=L) K = @. sqrt(grid.Krsq) # a 2D array with the total wavenumber @@ -132,7 +132,7 @@ fig # ## Setting initial conditions # Our initial condition is simply fluid at rest. -SingleLayerQG.set_q!(prob, ArrayType(dev)(zeros(grid.nx, grid.ny))) +SingleLayerQG.set_q!(prob, device_array(dev)(zeros(grid.nx, grid.ny))) # ## Diagnostics diff --git a/examples/singlelayerqg_decaying_topography.jl b/examples/singlelayerqg_decaying_topography.jl index 0a320ba8..389d4faf 100644 --- a/examples/singlelayerqg_decaying_topography.jl +++ b/examples/singlelayerqg_decaying_topography.jl @@ -88,7 +88,7 @@ fig # Our initial condition consist of a flow that has power only at wavenumbers with # ``6 < \frac{L}{2\pi} \sqrt{k_x^2 + k_y^2} < 12`` and initial energy ``E_0``. -# `ArrayType()` function returns the array type appropriate for the device, i.e., `Array` for +# `device_array()` function returns the array type appropriate for the device, i.e., `Array` for # `dev = CPU()` and `CuArray` for `dev = GPU()`. E₀ = 0.04 # energy of initial condition @@ -96,7 +96,7 @@ E₀ = 0.04 # energy of initial condition K = @. sqrt(grid.Krsq) # a 2D array with the total wavenumber Random.seed!(1234) -qih = ArrayType(dev)(randn(Complex{eltype(grid)}, size(sol))) +qih = device_array(dev)(randn(Complex{eltype(grid)}, size(sol))) @. qih = ifelse(K < 6 * 2π/L, 0, qih) @. qih = ifelse(K > 12 * 2π/L, 0, qih) qih *= sqrt(E₀ / SingleLayerQG.energy(qih, vars, params, grid)) # normalize qi to have energy E₀ diff --git a/examples/twodnavierstokes_stochasticforcing.jl b/examples/twodnavierstokes_stochasticforcing.jl index b0fc5abd..188127c2 100644 --- a/examples/twodnavierstokes_stochasticforcing.jl +++ b/examples/twodnavierstokes_stochasticforcing.jl @@ -55,7 +55,7 @@ forcing_wavenumber = 14.0 * 2π/L # the forcing wavenumber, `k_f`, for a spectr forcing_bandwidth = 1.5 * 2π/L # the width of the forcing spectrum, `δ_f` ε = 0.1 # energy input rate by the forcing -grid = TwoDGrid(dev, n, L) +grid = TwoDGrid(dev; nx=n, Lx=L) K = @. sqrt(grid.Krsq) # a 2D array with the total wavenumber @@ -125,7 +125,7 @@ fig # ## Setting initial conditions # Our initial condition is a fluid at rest. -TwoDNavierStokes.set_ζ!(prob, ArrayType(dev)(zeros(grid.nx, grid.ny))) +TwoDNavierStokes.set_ζ!(prob, device_array(dev)(zeros(grid.nx, grid.ny))) # ## Diagnostics diff --git a/examples/twodnavierstokes_stochasticforcing_budgets.jl b/examples/twodnavierstokes_stochasticforcing_budgets.jl index 7707f475..96845770 100644 --- a/examples/twodnavierstokes_stochasticforcing_budgets.jl +++ b/examples/twodnavierstokes_stochasticforcing_budgets.jl @@ -57,7 +57,7 @@ forcing_wavenumber = 14.0 * 2π/L # the forcing wavenumber, `k_f`, for a spectr forcing_bandwidth = 1.5 * 2π/L # the width of the forcing spectrum, `δ_f` ε = 0.1 # energy input rate by the forcing -grid = TwoDGrid(dev, n, L) +grid = TwoDGrid(dev; nx=n, Lx=L) K = @. sqrt(grid.Krsq) # a 2D array with the total wavenumber @@ -128,7 +128,7 @@ fig # ## Setting initial conditions # Our initial condition is a fluid at rest. -TwoDNavierStokes.set_ζ!(prob, ArrayType(dev)(zeros(grid.nx, grid.ny))) +TwoDNavierStokes.set_ζ!(prob, device_array(dev)(zeros(grid.nx, grid.ny))) # ## Diagnostics diff --git a/src/barotropicqgql.jl b/src/barotropicqgql.jl index b3e8963a..e2213d43 100644 --- a/src/barotropicqgql.jl +++ b/src/barotropicqgql.jl @@ -90,21 +90,21 @@ function Problem(dev::Device=CPU(); T = Float64) # the grid - grid = TwoDGrid(dev, nx, Lx, ny, Ly; aliased_fraction=aliased_fraction, T=T) + grid = TwoDGrid(dev; nx, Lx, ny, Ly, aliased_fraction=aliased_fraction, T) x, y = gridpoints(grid) # topographic PV eta === nothing && ( eta = zeros(dev, T, (nx, ny)) ) - params = !(typeof(eta)<:ArrayType(dev)) ? + params = !(typeof(eta)<:device_array(dev)) ? Params(grid, β, eta, μ, ν, nν, calcF) : Params(β, eta, rfft(eta), μ, ν, nν, calcF) - vars = calcF == nothingfunction ? DecayingVars(dev, grid) : stochastic ? StochasticForcedVars(dev, grid) : ForcedVars(dev, grid) + vars = calcF == nothingfunction ? DecayingVars(grid) : stochastic ? StochasticForcedVars(grid) : ForcedVars(grid) equation = BarotropicQGQL.Equation(params, grid) - FourierFlows.Problem(equation, stepper, dt, grid, vars, params, dev) + FourierFlows.Problem(equation, stepper, dt, grid, vars, params) end @@ -232,13 +232,14 @@ const ForcedVars = Vars{<:AbstractArray, <:AbstractArray, <:AbstractArray, Nothi const StochasticForcedVars = Vars{<:AbstractArray, <:AbstractArray, <:AbstractArray, <:AbstractArray} """ - DecayingVars(dev, grid) + DecayingVars(grid) -Return the vars for unforced two-dimensional quasi-linear barotropic QG problem on device `dev` -and with `grid`. +Return the vars for unforced two-dimensional quasi-linear barotropic QG problem on `grid`. """ -function DecayingVars(dev::Dev, grid::AbstractGrid) where Dev +function DecayingVars(grid::AbstractGrid) + Dev = typeof(grid.device) T = eltype(grid) + @devzeros Dev T (grid.nx, grid.ny) u v U uzeta vzeta zeta Zeta psi Psi @devzeros Dev Complex{T} (grid.nkr, grid.nl) N NZ uh vh Uh zetah Zetah psih Psih @@ -246,13 +247,14 @@ function DecayingVars(dev::Dev, grid::AbstractGrid) where Dev end """ - ForcedVars(dev, grid) + ForcedVars(grid) -Return the `vars` for forced two-dimensional quasi-linear barotropic QG problem on device -`dev` and with `grid`. +Return the `vars` for forced two-dimensional quasi-linear barotropic QG problem on `grid`. """ -function ForcedVars(dev::Dev, grid::AbstractGrid) where Dev +function ForcedVars(grid::AbstractGrid) + Dev = typeof(grid.device) T = eltype(grid) + @devzeros Dev T (grid.nx, grid.ny) u v U uzeta vzeta zeta Zeta psi Psi @devzeros Dev Complex{T} (grid.nkr, grid.nl) N NZ uh vh Uh zetah Zetah psih Psih Fh @@ -260,16 +262,18 @@ function ForcedVars(dev::Dev, grid::AbstractGrid) where Dev end """ - StochasticForcedVars(dev, grid) + StochasticForcedVars(grid) Return the `vars` for stochastically forced two-dimensional quasi-linear barotropic QG problem -on device `dev` and with `grid`. +on `grid`. """ -function StochasticForcedVars(dev::Dev, grid::AbstractGrid) where Dev +function StochasticForcedVars(grid::AbstractGrid) + Dev = typeof(grid.device) T = eltype(grid) + @devzeros Dev T (grid.nx, grid.ny) u v U uzeta vzeta zeta Zeta psi Psi @devzeros Dev Complex{T} (grid.nkr, grid.nl) N NZ uh vh Uh zetah Zetah psih Psih Fh prevsol - + return Vars(u, v, U, uzeta, vzeta, zeta, Zeta, psi, Psi, N, NZ, uh, vh, Uh, zetah, Zetah, psih, Psih, Fh, prevsol) end @@ -335,11 +339,11 @@ N = - \\widehat{𝖩(ψ, ζ + η)}^{\\mathrm{QL}} + F̂ . """ function calcN!(N, sol, t, clock, vars, params, grid) dealias!(sol, grid) - + calcN_advection!(N, sol, t, clock, vars, params, grid) - + addforcing!(N, sol, t, clock, vars, params, grid) - + return nothing end diff --git a/src/multilayerqg.jl b/src/multilayerqg.jl index 8475d8b8..f4b786ed 100644 --- a/src/multilayerqg.jl +++ b/src/multilayerqg.jl @@ -121,15 +121,15 @@ function Problem(nlayers::Int, # number of fluid layers # topographic PV eta === nothing && (eta = zeros(dev, T, (nx, ny))) - grid = TwoDGrid(dev, nx, Lx, ny, Ly; aliased_fraction=aliased_fraction, T=T) + grid = TwoDGrid(dev; nx, Lx, ny, Ly, aliased_fraction, T) - params = Params(nlayers, g, f₀, β, ρ, H, U, eta, μ, ν, nν, grid, calcFq=calcFq, dev=dev) + params = Params(nlayers, g, f₀, β, ρ, H, U, eta, μ, ν, nν, grid, calcFq=calcFq) - vars = calcFq == nothingfunction ? DecayingVars(dev, grid, params) : (stochastic ? StochasticForcedVars(dev, grid, params) : ForcedVars(dev, grid, params)) + vars = calcFq == nothingfunction ? DecayingVars(grid, params) : (stochastic ? StochasticForcedVars(grid, params) : ForcedVars(grid, params)) - equation = linear ? LinearEquation(dev, params, grid) : Equation(dev, params, grid) + equation = linear ? LinearEquation(params, grid) : Equation(params, grid) - FourierFlows.Problem(equation, stepper, dt, grid, vars, params, dev) + FourierFlows.Problem(equation, stepper, dt, grid, vars, params) end """ @@ -281,14 +281,15 @@ end function convert_U_to_U3D(dev, nlayers, grid, U::Number) T = eltype(grid) - A = ArrayType(dev) + A = device_array(dev) U_3D = reshape(repeat([T(U)], outer=(grid.ny, 1)), (1, grid.ny, nlayers)) return A(U_3D) end -function Params(nlayers, g, f₀, β, ρ, H, U, eta, μ, ν, nν, grid; calcFq=nothingfunction, effort=FFTW.MEASURE, dev::Device=CPU()) where TU +function Params(nlayers, g, f₀, β, ρ, H, U, eta, μ, ν, nν, grid; calcFq=nothingfunction, effort=FFTW.MEASURE) where TU + dev = grid.device T = eltype(grid) - A = ArrayType(dev) + A = device_array(dev) ny, nx = grid.ny , grid.nx nkr, nl = grid.nkr, grid.nl @@ -358,7 +359,7 @@ numberoflayers(::TwoLayerParams) = 2 # --------- """ - hyperviscosity(dev, params, grid) + hyperviscosity(params, grid) Return the linear operator `L` that corresponds to (hyper)-viscosity of order ``n_ν`` with coefficient ``ν`` for ``n`` fluid layers. @@ -366,9 +367,11 @@ coefficient ``ν`` for ``n`` fluid layers. L_j = - ν |𝐤|^{2 n_ν}, \\ j = 1, ...,n . ``` """ -function hyperviscosity(dev, params, grid) +function hyperviscosity(params, grid) + dev = grid.device T = eltype(grid) - L = ArrayType(dev){T}(undef, (grid.nkr, grid.nl, numberoflayers(params))) + + L = device_array(dev){T}(undef, (grid.nkr, grid.nl, numberoflayers(params))) @. L = - params.ν * grid.Krsq^params.nν @views @. L[1, 1, :] = 0 @@ -376,31 +379,31 @@ function hyperviscosity(dev, params, grid) end """ - LinearEquation(dev, params, grid) + LinearEquation(params, grid) Return the equation for a multi-layer quasi-geostrophic problem with `params` and `grid`. The linear opeartor ``L`` includes only (hyper)-viscosity and is computed via -`hyperviscosity(dev, params, grid)`. +`hyperviscosity(params, grid)`. The nonlinear term is computed via function `calcNlinear!`. """ -function LinearEquation(dev, params, grid) - L = hyperviscosity(dev, params, grid) +function LinearEquation(params, grid) + L = hyperviscosity(params, grid) return FourierFlows.Equation(L, calcNlinear!, grid) end """ - Equation(dev, params, grid) + Equation(params, grid) Return the equation for a multi-layer quasi-geostrophic problem with `params` and `grid`. The linear opeartor ``L`` includes only (hyper)-viscosity and is computed via -`hyperviscosity(dev, params, grid)`. +`hyperviscosity(params, grid)`. The nonlinear term is computed via function `calcN!`. """ -function Equation(dev, params, grid) - L = hyperviscosity(dev, params, grid) +function Equation(params, grid) + L = hyperviscosity(params, grid) return FourierFlows.Equation(L, calcN!, grid) end @@ -445,11 +448,12 @@ const ForcedVars = Vars{<:AbstractArray, <:AbstractArray, <:AbstractArray, Nothi const StochasticForcedVars = Vars{<:AbstractArray, <:AbstractArray, <:AbstractArray, <:AbstractArray} """ - DecayingVars(dev, grid, params) + DecayingVars(grid, params) Return the variables for an unforced multi-layer QG problem with `grid` and `params`. """ -function DecayingVars(dev::Dev, grid, params) where Dev +function DecayingVars(grid, params) + Dev = typeof(grid.device) T = eltype(grid) nlayers = numberoflayers(params) @@ -460,11 +464,12 @@ function DecayingVars(dev::Dev, grid, params) where Dev end """ - ForcedVars(dev, grid, params) + ForcedVars(grid, params) Return the variables for a forced multi-layer QG problem with `grid` and `params`. """ -function ForcedVars(dev::Dev, grid, params) where Dev +function ForcedVars(grid, params) + Dev = typeof(grid.device) T = eltype(grid) nlayers = numberoflayers(params) @@ -475,11 +480,12 @@ function ForcedVars(dev::Dev, grid, params) where Dev end """ - StochasticForcedVars(dev, rid, params) + StochasticForcedVars(rid, params) Return the variables for a forced multi-layer QG problem with `grid` and `params`. """ -function StochasticForcedVars(dev::Dev, grid, params) where Dev +function StochasticForcedVars(grid, params) + Dev = typeof(grid.device) T = eltype(grid) nlayers = numberoflayers(params) diff --git a/test/test_barotropicqgql.jl b/test/test_barotropicqgql.jl index 4f4362aa..a0133133 100644 --- a/test/test_barotropicqgql.jl +++ b/test/test_barotropicqgql.jl @@ -61,7 +61,7 @@ function test_bqgql_stochasticforcingbudgets(dev::Device=CPU(); n=256, dt=0.01, kf, dkf = 12.0, 2.0 ε = 0.1 - Kr = CUDA.@allowscalar ArrayType(dev)([ grid.kr[i] for i=1:grid.nkr, j=1:grid.nl]) + Kr = CUDA.@allowscalar device_array(dev)([ grid.kr[i] for i=1:grid.nkr, j=1:grid.nl]) forcing_spectrum = zeros(dev, T, (grid.nkr, grid.nl)) @. forcing_spectrum = exp(-(sqrt(grid.Krsq) - kf)^2 / (2 * dkf^2)) @@ -74,7 +74,7 @@ function test_bqgql_stochasticforcingbudgets(dev::Device=CPU(); n=256, dt=0.01, Random.seed!(1234) function calcF!(F, sol, t, clock, vars, params, grid) - eta = ArrayType(dev)(exp.(2π * im * rand(T, size(sol))) / sqrt(clock.dt)) + eta = device_array(dev)(exp.(2π * im * rand(T, size(sol))) / sqrt(clock.dt)) CUDA.@allowscalar eta[1, 1] = 0 @. F = eta * sqrt(forcing_spectrum) return nothing @@ -230,7 +230,7 @@ end function test_bqgql_problemtype(dev, T) prob = BarotropicQGQL.Problem(dev; T=T) - A = ArrayType(dev) + 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}) end diff --git a/test/test_multilayerqg.jl b/test/test_multilayerqg.jl index 5f973ddd..3549859f 100644 --- a/test/test_multilayerqg.jl +++ b/test/test_multilayerqg.jl @@ -165,7 +165,7 @@ that a solution to the problem forced by this Ff is then qf. function test_mqg_nonlinearadvection(dt, stepper, dev::Device=CPU(); n=128, L=2π, nlayers=2, μ=0.0, ν=0.0, nν=1) - A = ArrayType(dev) + A = device_array(dev) tf = 0.5 nt = round(Int, tf/dt) @@ -256,7 +256,7 @@ is unstable.) function test_mqg_linearadvection(dt, stepper, dev::Device=CPU(); n=128, L=2π, nlayers=2, μ=0.0, ν=0.0, nν=1) - A = ArrayType(dev) + A = device_array(dev) tf = 0.5 nt = round(Int, tf/dt) @@ -550,7 +550,7 @@ function test_mqg_problemtype(dev, T) prob1 = MultiLayerQG.Problem(1, dev; T) prob2 = MultiLayerQG.Problem(2, dev; T) - A = ArrayType(dev) + A = device_array(dev) return typeof(prob1.sol)<:A{Complex{T}, 3} && typeof(prob1.grid.Lx)==T && diff --git a/test/test_singlelayerqg.jl b/test/test_singlelayerqg.jl index 47825381..642b03c7 100644 --- a/test/test_singlelayerqg.jl +++ b/test/test_singlelayerqg.jl @@ -73,7 +73,7 @@ function test_1layerqg_stochasticforcing_energybudget(dev::Device=CPU(); n=256, Random.seed!(1234) function calcF!(Fh, sol, t, clock, vars, params, grid) - eta = ArrayType(dev)(exp.(2π * im * rand(T, size(sol))) / sqrt(clock.dt)) + eta = device_array(dev)(exp.(2π * im * rand(T, size(sol))) / sqrt(clock.dt)) CUDA.@allowscalar eta[1, 1] = 0 @. Fh = eta * sqrt(forcing_spectrum) @@ -179,7 +179,7 @@ function test_1layerqg_stochasticforcing_enstrophybudget(dev::Device=CPU(); n=25 Random.seed!(1234) function calcF!(Fh, sol, t, clock, vars, params, grid) - eta = ArrayType(dev)(exp.(2π * im * rand(T, size(sol))) / sqrt(clock.dt)) + eta = device_array(dev)(exp.(2π * im * rand(T, size(sol))) / sqrt(clock.dt)) CUDA.@allowscalar eta[1, 1] = 0 @. Fh = eta * sqrt(forcing_spectrum) @@ -373,7 +373,7 @@ Tests 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) - A = ArrayType(dev) + 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}) end diff --git a/test/test_surfaceqg.jl b/test/test_surfaceqg.jl index a9a74542..6cf79ad1 100644 --- a/test/test_surfaceqg.jl +++ b/test/test_surfaceqg.jl @@ -1,7 +1,7 @@ function test_sqg_problemtype(dev, T) prob = SurfaceQG.Problem(dev; T=T) - A = ArrayType(dev) + A = device_array(dev) (typeof(prob.sol)<:A{Complex{T},2} && typeof(prob.grid.Lx)==T && eltype(prob.grid.x)==T && typeof(prob.vars.u)<:A{T,2}) end @@ -115,7 +115,7 @@ function test_sqg_noforcing(dev::Device=CPU()) SurfaceQG.addforcing!(prob_unforced.timestepper.N, prob_unforced.sol, prob_unforced.clock.t, prob_unforced.clock, prob_unforced.vars, prob_unforced.params, prob_unforced.grid) function calcF!(Fh, sol, t, clock, vars, params, grid) - Fh .= 2 * ArrayType(dev)(ones(size(sol))) + Fh .= 2 * device_array(dev)(ones(size(sol))) return nothing end @@ -123,7 +123,7 @@ function test_sqg_noforcing(dev::Device=CPU()) SurfaceQG.addforcing!(prob_forced.timestepper.N, prob_forced.sol, prob_forced.clock.t, prob_forced.clock, prob_forced.vars, prob_forced.params, prob_forced.grid) - return prob_unforced.timestepper.N == Complex.(ArrayType(dev)(zeros(size(prob_unforced.sol)))) && prob_forced.timestepper.N == Complex.(2*ArrayType(dev)(ones(size(prob_unforced.sol)))) + return prob_unforced.timestepper.N == Complex.(device_array(dev)(zeros(size(prob_unforced.sol)))) && prob_forced.timestepper.N == Complex.(2*device_array(dev)(ones(size(prob_unforced.sol)))) end function test_sqg_deterministicforcing_buoyancy_variance_budget(dev::Device=CPU(); n=256, dt=0.01, L=2π, ν=1e-7, nν=2, tf=10.0) @@ -173,9 +173,9 @@ function test_sqg_stochasticforcing_buoyancy_variance_budget(dev::Device=CPU(); grid = TwoDGrid(dev, n, L) x, y = gridpoints(grid) - Kr = ArrayType(dev)([CUDA.@allowscalar grid.kr[i] for i=1:grid.nkr, j=1:grid.nl]) + Kr = device_array(dev)([CUDA.@allowscalar grid.kr[i] for i=1:grid.nkr, j=1:grid.nl]) - forcing_spectrum = ArrayType(dev)(zero(grid.Krsq)) + forcing_spectrum = device_array(dev)(zero(grid.Krsq)) @. forcing_spectrum = exp(-(sqrt(grid.Krsq) - kf)^2 / (2 * dkf^2)) @. forcing_spectrum = ifelse(grid.Krsq < 2^2, 0, forcing_spectrum) @. forcing_spectrum = ifelse(grid.Krsq > 20^2, 0, forcing_spectrum) @@ -186,7 +186,7 @@ function test_sqg_stochasticforcing_buoyancy_variance_budget(dev::Device=CPU(); Random.seed!(1234) function calcF!(Fh, sol, t, clock, vars, params, grid) - eta = ArrayType(dev)(exp.(2π * im * rand(Float64, size(sol))) / sqrt(clock.dt)) + eta = device_array(dev)(exp.(2π * im * rand(Float64, size(sol))) / sqrt(clock.dt)) CUDA.@allowscalar eta[1, 1] = 0.0 @. Fh = eta * sqrt(forcing_spectrum) nothing diff --git a/test/test_twodnavierstokes.jl b/test/test_twodnavierstokes.jl index 12a2b07e..050cd3ee 100644 --- a/test/test_twodnavierstokes.jl +++ b/test/test_twodnavierstokes.jl @@ -32,9 +32,9 @@ function test_twodnavierstokes_stochasticforcing_energybudget(dev::Device=CPU(); grid = TwoDGrid(dev, n, L) x, y = gridpoints(grid) - Kr = ArrayType(dev)([CUDA.@allowscalar grid.kr[i] for i=1:grid.nkr, j=1:grid.nl]) + Kr = device_array(dev)([CUDA.@allowscalar grid.kr[i] for i=1:grid.nkr, j=1:grid.nl]) - forcing_spectrum = ArrayType(dev)(zero(grid.Krsq)) + forcing_spectrum = device_array(dev)(zero(grid.Krsq)) @. forcing_spectrum = exp(-(sqrt(grid.Krsq) - kf)^2 / (2 * dkf^2)) @. forcing_spectrum = ifelse(grid.Krsq < 2^2, 0, forcing_spectrum) @. forcing_spectrum = ifelse(grid.Krsq > 20^2, 0, forcing_spectrum) @@ -45,7 +45,7 @@ function test_twodnavierstokes_stochasticforcing_energybudget(dev::Device=CPU(); Random.seed!(1234) function calcF!(Fh, sol, t, clock, vars, params, grid) - eta = ArrayType(dev)(exp.(2π * im * rand(Float64, size(sol))) / sqrt(clock.dt)) + eta = device_array(dev)(exp.(2π * im * rand(Float64, size(sol))) / sqrt(clock.dt)) CUDA.@allowscalar eta[1, 1] = 0.0 @. Fh = eta * sqrt(forcing_spectrum) @@ -84,9 +84,9 @@ function test_twodnavierstokes_stochasticforcing_enstrophybudget(dev::Device=CPU grid = TwoDGrid(dev, n, L) x, y = gridpoints(grid) - Kr = ArrayType(dev)([CUDA.@allowscalar grid.kr[i] for i=1:grid.nkr, j=1:grid.nl]) + Kr = device_array(dev)([CUDA.@allowscalar grid.kr[i] for i=1:grid.nkr, j=1:grid.nl]) - forcing_spectrum = ArrayType(dev)(zero(grid.Krsq)) + forcing_spectrum = device_array(dev)(zero(grid.Krsq)) @. forcing_spectrum = exp(-(sqrt(grid.Krsq) - kf)^2 / (2 * dkf^2)) @. forcing_spectrum = ifelse(grid.Krsq < 2^2, 0, forcing_spectrum) @. forcing_spectrum = ifelse(grid.Krsq > 20^2, 0, forcing_spectrum) @@ -97,7 +97,7 @@ function test_twodnavierstokes_stochasticforcing_enstrophybudget(dev::Device=CPU Random.seed!(1234) function calcF!(Fh, sol, t, cl, v, p, g) - eta = ArrayType(dev)(exp.(2π * im * rand(Float64, size(sol))) / sqrt(cl.dt)) + eta = device_array(dev)(exp.(2π * im * rand(Float64, size(sol))) / sqrt(cl.dt)) CUDA.@allowscalar eta[1, 1] = 0.0 @. Fh = eta * sqrt(forcing_spectrum) @@ -291,7 +291,7 @@ end function test_twodnavierstokes_problemtype(dev, T) prob = TwoDNavierStokes.Problem(dev; T=T) - A = ArrayType(dev) + A = device_array(dev) (typeof(prob.sol)<:A{Complex{T},2} && typeof(prob.grid.Lx)==T && eltype(prob.grid.x)==T && typeof(prob.vars.u)<:A{T,2}) end From 61258478283bcf0bc1de14215a36f8dc5d3f87f9 Mon Sep 17 00:00:00 2001 From: "Navid C. Constantinou" Date: Fri, 2 Sep 2022 13:48:42 +1000 Subject: [PATCH 2/7] update to latest grid refactor --- test/test_barotropicqgql.jl | 6 +++--- test/test_multilayerqg.jl | 4 ++-- test/test_singlelayerqg.jl | 10 +++++----- test/test_surfaceqg.jl | 6 +++--- test/test_twodnavierstokes.jl | 10 +++++----- test/test_utils.jl | 2 +- 6 files changed, 19 insertions(+), 19 deletions(-) diff --git a/test/test_barotropicqgql.jl b/test/test_barotropicqgql.jl index a0133133..9a5b4bfd 100644 --- a/test/test_barotropicqgql.jl +++ b/test/test_barotropicqgql.jl @@ -54,7 +54,7 @@ function test_bqgql_stochasticforcingbudgets(dev::Device=CPU(); n=256, dt=0.01, dt, tf = 0.005, 0.1/μ nt = round(Int, tf/dt) - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) # Forcing @@ -113,7 +113,7 @@ function test_bqgql_deterministicforcingbudgets(dev::Device=CPU(); n=256, dt=0.0 dt, tf = 0.005, 0.1/μ nt = round(Int, tf/dt) - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) k₀, l₀ = 2π/grid.Lx, 2π/grid.Ly @@ -166,7 +166,7 @@ function test_bqgql_advection(dt, stepper, dev::Device=CPU(); n=128, L=2π, ν=1 tf = 1.0 nt = round(Int, tf/dt) - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) ψf = @. cos(3y) + sin(2x)*cos(2y) + 2sin(x)*cos(3y) diff --git a/test/test_multilayerqg.jl b/test/test_multilayerqg.jl index 3549859f..12091091 100644 --- a/test/test_multilayerqg.jl +++ b/test/test_multilayerqg.jl @@ -73,7 +73,7 @@ q1 and q2. Similarly, that streamfunctionfrompv gives ψ1 and ψ2 from q1 and q2 """ function test_pvtofromstreamfunction_2layer(dev::Device=CPU()) n, L = 128, 2π - gr = TwoDGrid(dev, n, L) + gr = TwoDGrid(dev; nx=n, Lx=L) nlayers = 2 # these choice of parameters give the f₀, g = 1, 1 # desired PV-streamfunction relations @@ -116,7 +116,7 @@ q1, q2, and q3. """ function test_pvtofromstreamfunction_3layer(dev::Device=CPU()) n, L = 128, 2π - gr = TwoDGrid(dev, n, L) + gr = TwoDGrid(dev; nx=n, Lx=L) nlayers = 3 # these choice of parameters give the f₀, g = 1, 1 # desired PV-streamfunction relations diff --git a/test/test_singlelayerqg.jl b/test/test_singlelayerqg.jl index 642b03c7..cba25c62 100644 --- a/test/test_singlelayerqg.jl +++ b/test/test_singlelayerqg.jl @@ -55,7 +55,7 @@ function test_1layerqg_stochasticforcing_energybudget(dev::Device=CPU(); n=256, dt, tf = 0.005, 0.1/μ nt = round(Int, tf/dt) - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) K = @. sqrt(grid.Krsq) # a 2D array with the total wavenumber @@ -114,7 +114,7 @@ function test_1layerqg_deterministicforcing_energybudget(dev::Device=CPU(); n=25 dt, tf = 0.005, 0.1/μ nt = round(Int, tf/dt) - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) k₀, l₀ = 2π / grid.Lx, 2π/grid.Ly @@ -161,7 +161,7 @@ function test_1layerqg_stochasticforcing_enstrophybudget(dev::Device=CPU(); n=25 dt, tf = 0.005, 0.1/μ nt = round(Int, tf/dt) - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) K = @. sqrt(grid.Krsq) # a 2D array with the total wavenumber @@ -220,7 +220,7 @@ function test_1layerqg_deterministicforcing_enstrophybudget(dev::Device=CPU(); n dt, tf = 0.005, 0.1/μ nt = round(Int, tf/dt) - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) k₀, l₀ = 2π/grid.Lx, 2π/grid.Ly @@ -272,7 +272,7 @@ function test_1layerqg_advection(dt, stepper, dev::Device=CPU(); n=128, L=2π, tf = 1.0 nt = round(Int, tf/dt) - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) ψf = @. sin(2x) * cos(2y) + 2sin(x) * cos(3y) diff --git a/test/test_surfaceqg.jl b/test/test_surfaceqg.jl index 6cf79ad1..7de12de1 100644 --- a/test/test_surfaceqg.jl +++ b/test/test_surfaceqg.jl @@ -35,7 +35,7 @@ function test_sqg_advection(dt, stepper, dev::Device=CPU(); n=128, L=2π, ν=1e- tf = 0.5 # SQG piles up energy at small scales so running for t ⪆ 0.5 brings instability nt = round(Int, tf/dt) - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) ψf = @. sin(2x)*cos(2y) + 2sin(x)*cos(3y) @@ -132,7 +132,7 @@ function test_sqg_deterministicforcing_buoyancy_variance_budget(dev::Device=CPU( dt, tf = 0.005, 10 nt = round(Int, tf/dt) - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) # buoyancy forcing = 0.01cos(4x)cos(5y)cos(2t) @@ -170,7 +170,7 @@ function test_sqg_stochasticforcing_buoyancy_variance_budget(dev::Device=CPU(); kf, dkf = 12.0, 2.0 εᵇ = 0.01 - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) Kr = device_array(dev)([CUDA.@allowscalar grid.kr[i] for i=1:grid.nkr, j=1:grid.nl]) diff --git a/test/test_twodnavierstokes.jl b/test/test_twodnavierstokes.jl index 050cd3ee..209b2571 100644 --- a/test/test_twodnavierstokes.jl +++ b/test/test_twodnavierstokes.jl @@ -29,7 +29,7 @@ function test_twodnavierstokes_stochasticforcing_energybudget(dev::Device=CPU(); kf, dkf = 12.0, 2.0 ε = 0.1 - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) Kr = device_array(dev)([CUDA.@allowscalar grid.kr[i] for i=1:grid.nkr, j=1:grid.nl]) @@ -81,7 +81,7 @@ function test_twodnavierstokes_stochasticforcing_enstrophybudget(dev::Device=CPU kf, dkf = 12.0, 2.0 εᶻ = 0.1 - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) Kr = device_array(dev)([CUDA.@allowscalar grid.kr[i] for i=1:grid.nkr, j=1:grid.nl]) @@ -133,7 +133,7 @@ function test_twodnavierstokes_deterministicforcing_energybudget(dev::Device=CPU dt, tf = 0.005, 0.1/μ nt = round(Int, tf/dt) - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) # Forcing = 0.01cos(4x)cos(5y)cos(2t) @@ -174,7 +174,7 @@ function test_twodnavierstokes_deterministicforcing_enstrophybudget(dev::Device= dt, tf = 0.005, 0.1/μ nt = round(Int, tf/dt) - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) # Forcing = 0.01cos(4x)cos(5y)cos(2t) @@ -226,7 +226,7 @@ function test_twodnavierstokes_advection(dt, stepper, dev::Device=CPU(); n=128, tf = 1.0 nt = round(Int, tf/dt) - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) x, y = gridpoints(grid) ψf = @. sin(2x)*cos(2y) + 2sin(x)*cos(3y) diff --git a/test/test_utils.jl b/test/test_utils.jl index 54f1d2bc..7724e1e0 100644 --- a/test/test_utils.jl +++ b/test/test_utils.jl @@ -3,7 +3,7 @@ Test the peakedisotropicspectrum function. """ function testpeakedisotropicspectrum(dev::Device=CPU()) n, L = 128, 2π - grid = TwoDGrid(dev, n, L) + grid = TwoDGrid(dev; nx=n, Lx=L) k0, E0 = 6, 0.5 qi = peakedisotropicspectrum(grid, k0, E0; allones=true) ρ, qhρ = FourierFlows.radialspectrum(rfft(qi) .* grid.invKrsq, grid) From ad0dbc5b6e6704cba0e3e64d361ffb6ae3f4d717 Mon Sep 17 00:00:00 2001 From: "Navid C. Constantinou" Date: Fri, 2 Sep 2022 13:51:47 +1000 Subject: [PATCH 3/7] update to latest grid refactor --- src/singlelayerqg.jl | 28 ++++++++++++++++------------ 1 file changed, 16 insertions(+), 12 deletions(-) diff --git a/src/singlelayerqg.jl b/src/singlelayerqg.jl index 60cfd34a..4d1caab6 100644 --- a/src/singlelayerqg.jl +++ b/src/singlelayerqg.jl @@ -93,7 +93,7 @@ function Problem(dev::Device=CPU(); T = Float64) # the grid - grid = TwoDGrid(dev, nx, Lx, ny, Ly; aliased_fraction=aliased_fraction, T=T) + grid = TwoDGrid(dev; nx, Lx, ny, Ly, aliased_fraction, T) x, y = gridpoints(grid) # topographic PV @@ -101,11 +101,11 @@ function Problem(dev::Device=CPU(); params = deformation_radius == Inf ? BarotropicQGParams(grid, T(β), eta, T(μ), T(ν), nν, calcF) : EquivalentBarotropicQGParams(grid, T(β), T(deformation_radius), eta, T(μ), T(ν), nν, calcF) - vars = calcF == nothingfunction ? DecayingVars(dev, grid) : (stochastic ? StochasticForcedVars(dev, grid) : ForcedVars(dev, grid)) + vars = calcF == nothingfunction ? DecayingVars(grid) : (stochastic ? StochasticForcedVars(grid) : ForcedVars(grid)) equation = Equation(params, grid) - return FourierFlows.Problem(equation, stepper, dt, grid, vars, params, dev) + return FourierFlows.Problem(equation, stepper, dt, grid, vars, params) end @@ -250,11 +250,12 @@ const ForcedVars = Vars{<:AbstractArray, <:AbstractArray, <:AbstractArray, Nothi const StochasticForcedVars = Vars{<:AbstractArray, <:AbstractArray, <:AbstractArray, <:AbstractArray} """ - DecayingVars(dev, grid) + DecayingVars(grid) -Return the `vars` for unforced single-layer QG problem on device `dev` and with `grid` +Return the `vars` for unforced single-layer QG problem on `grid`. """ -function DecayingVars(dev::Dev, grid::AbstractGrid) where Dev +function DecayingVars(grid::AbstractGrid) + Dev = typeof(grid.device) T = eltype(grid) @devzeros Dev T (grid.nx, grid.ny) q u v ψ @@ -264,11 +265,12 @@ function DecayingVars(dev::Dev, grid::AbstractGrid) where Dev end """ - ForcedVars(dev, grid) + ForcedVars(grid) -Return the `vars` for forced single-layer QG problem on device dev and with `grid`. +Return the `vars` for forced single-layer QG problem on `grid`. """ -function ForcedVars(dev::Dev, grid::AbstractGrid) where Dev +function ForcedVars(grid::AbstractGrid) + Dev = typeof(grid.device) T = eltype(grid) @devzeros Dev T (grid.nx, grid.ny) q u v ψ @@ -278,12 +280,14 @@ function ForcedVars(dev::Dev, grid::AbstractGrid) where Dev end """ - StochasticForcedVars(dev, grid) + StochasticForcedVars(grid) -Return the vars for stochastically forced barotropic QG problem on device dev and with `grid`. +Return the vars for stochastically forced barotropic QG problem on `grid`. """ -function StochasticForcedVars(dev::Dev, grid::AbstractGrid) where Dev +function StochasticForcedVars(grid::AbstractGrid) + Dev = typeof(grid.device) T = eltype(grid) + @devzeros Dev T (grid.nx, grid.ny) q u v ψ @devzeros Dev Complex{T} (grid.nkr, grid.nl) qh uh vh ψh Fh prevsol From 1e2d1a9b33aeab97b67d941f6993c18abf44abb4 Mon Sep 17 00:00:00 2001 From: "Navid C. Constantinou" Date: Fri, 2 Sep 2022 13:54:33 +1000 Subject: [PATCH 4/7] update to latest grid refactor --- src/surfaceqg.jl | 29 +++++++++++++++++------------ src/twodnavierstokes.jl | 29 +++++++++++++++-------------- 2 files changed, 32 insertions(+), 26 deletions(-) diff --git a/src/surfaceqg.jl b/src/surfaceqg.jl index 5b7bb77e..e68e7439 100644 --- a/src/surfaceqg.jl +++ b/src/surfaceqg.jl @@ -75,11 +75,11 @@ function Problem(dev::Device=CPU(); aliased_fraction = 1/3, T = Float64) - grid = TwoDGrid(dev, nx, Lx, ny, Ly; aliased_fraction=aliased_fraction, T=T) + grid = TwoDGrid(dev; nx, Lx, ny, Ly, aliased_fraction, T) params = Params{T}(ν, nν, calcF) - vars = calcF == nothingfunction ? DecayingVars(dev, grid) : (stochastic ? StochasticForcedVars(dev, grid) : ForcedVars(dev, grid)) + vars = calcF == nothingfunction ? DecayingVars(grid) : (stochastic ? StochasticForcedVars(grid) : ForcedVars(grid)) equation = Equation(params, grid) @@ -173,12 +173,14 @@ const ForcedVars = Vars{<:AbstractArray, <:AbstractArray, <:AbstractArray, Nothi const StochasticForcedVars = Vars{<:AbstractArray, <:AbstractArray, <:AbstractArray, <:AbstractArray} """ - DecayingVars(dev, grid) + DecayingVars(grid) -Return the `vars` for unforced surface QG dynamics on device `dev` and with `grid`. +Return the `vars` for unforced surface QG dynamics on `grid`. """ -function DecayingVars(::Dev, grid::AbstractGrid) where Dev +function DecayingVars(grid::AbstractGrid) + Dev = typeof(grid.device) T = eltype(grid) + @devzeros Dev T (grid.nx, grid.ny) b u v @devzeros Dev Complex{T} (grid.nkr, grid.nl) bh uh vh @@ -186,12 +188,14 @@ function DecayingVars(::Dev, grid::AbstractGrid) where Dev end """ - ForcedVars(dev, grid) + ForcedVars(grid) -Return the vars for forced surface QG dynamics on device `dev` and with `grid`. +Return the vars for forced surface QG dynamics on `grid`. """ -function ForcedVars(dev::Dev, grid) where Dev +function ForcedVars(grid) + Dev = typeof(grid.device) T = eltype(grid) + @devzeros Dev T (grid.nx, grid.ny) b u v @devzeros Dev Complex{T} (grid.nkr, grid.nl) bh uh vh Fh @@ -199,13 +203,14 @@ function ForcedVars(dev::Dev, grid) where Dev end """ - StochasticForcedVars(dev, grid) + StochasticForcedVars(grid) -Return the `vars` for stochastically forced surface QG dynamics on device `dev` and with `grid`. +Return the `vars` for stochastically forced surface QG dynamics on `grid`. """ -function StochasticForcedVars(dev::Dev, grid) where Dev +function StochasticForcedVars(grid) + Dev = typeof(grid.device) T = eltype(grid) - + @devzeros Dev T (grid.nx, grid.ny) b u v @devzeros Dev Complex{T} (grid.nkr, grid.nl) bh uh vh Fh prevsol diff --git a/src/twodnavierstokes.jl b/src/twodnavierstokes.jl index ff80ed92..bf9dbfd4 100644 --- a/src/twodnavierstokes.jl +++ b/src/twodnavierstokes.jl @@ -83,11 +83,11 @@ function Problem(dev::Device=CPU(); aliased_fraction = 1/3, T = Float64) - grid = TwoDGrid(dev, nx, Lx, ny, Ly; aliased_fraction=aliased_fraction, T) + grid = TwoDGrid(dev; nx, Lx, ny, Ly, aliased_fraction, T) params = Params(T(ν), nν, T(μ), nμ, calcF) - vars = calcF == nothingfunction ? DecayingVars(dev, grid) : (stochastic ? StochasticForcedVars(dev, grid) : ForcedVars(dev, grid)) + vars = calcF == nothingfunction ? DecayingVars(grid) : (stochastic ? StochasticForcedVars(grid) : ForcedVars(grid)) equation = Equation(params, grid) @@ -188,12 +188,12 @@ const StochasticForcedVars = Vars{<:AbstractArray, <:AbstractArray, <:AbstractAr """ DecayingVars(dev, grid) -Return the variables `vars` for unforced two-dimensional Navier-Stokes problem on device `dev` and -with `grid`. +Return the variables `vars` for unforced two-dimensional Navier-Stokes problem on `grid`. """ -function DecayingVars(::Dev, grid::AbstractGrid) where Dev +function DecayingVars(grid::AbstractGrid) + Dev = typeof(grid.device) T = eltype(grid) - + @devzeros Dev T (grid.nx, grid.ny) ζ u v @devzeros Dev Complex{T} (grid.nkr, grid.nl) ζh uh vh @@ -201,11 +201,12 @@ function DecayingVars(::Dev, grid::AbstractGrid) where Dev end """ - ForcedVars(dev, grid) + ForcedVars(grid) -Return the variables `vars` for forced two-dimensional Navier-Stokes on device `dev` and with `grid`. +Return the variables `vars` for forced two-dimensional Navier-Stokes on `grid`. """ -function ForcedVars(dev::Dev, grid::AbstractGrid) where Dev +function ForcedVars(grid::AbstractGrid) + Dev = typeof(grid.device) T = eltype(grid) @devzeros Dev T (grid.nx, grid.ny) ζ u v @@ -215,14 +216,14 @@ function ForcedVars(dev::Dev, grid::AbstractGrid) where Dev end """ - StochasticForcedVars(dev, grid) + StochasticForcedVars(grid) -Return the variables `vars` for stochastically forced two-dimensional Navier-Stokes on device `dev` and -with `grid`. +Return the variables `vars` for stochastically forced two-dimensional Navier-Stokes on device `grid`. """ -function StochasticForcedVars(dev::Dev, grid::AbstractGrid) where Dev +function StochasticForcedVars(grid::AbstractGrid) + Dev = typeof(grid.device) T = eltype(grid) - + @devzeros Dev T (grid.nx, grid.ny) ζ u v @devzeros Dev Complex{T} (grid.nkr, grid.nl) ζh uh vh Fh prevsol From f1f9db8b97940e628dba45babc3adc8d86e7e697 Mon Sep 17 00:00:00 2001 From: "Navid C. Constantinou" Date: Fri, 2 Sep 2022 14:00:59 +1000 Subject: [PATCH 5/7] update to latest grid refactor --- test/test_barotropicqgql.jl | 2 +- test/test_multilayerqg.jl | 16 ++++++++-------- test/test_singlelayerqg.jl | 4 ++-- test/test_surfaceqg.jl | 2 +- test/test_utils.jl | 6 +++--- 5 files changed, 15 insertions(+), 15 deletions(-) diff --git a/test/test_barotropicqgql.jl b/test/test_barotropicqgql.jl index 9a5b4bfd..f15e981e 100644 --- a/test/test_barotropicqgql.jl +++ b/test/test_barotropicqgql.jl @@ -204,7 +204,7 @@ Tests the energy and enstrophy function for a BarotropicQGQL problem. function test_bqgql_energyenstrophy(dev::Device=CPU()) nx, Lx = 64, 2π ny, Ly = 64, 3π - grid = TwoDGrid(dev, nx, Lx, ny, Ly) + grid = TwoDGrid(dev; nx, Lx, ny, Ly) k₀, l₀ = 2π/Lx, 2π/Ly # fundamental wavenumbers x, y = gridpoints(grid) diff --git a/test/test_multilayerqg.jl b/test/test_multilayerqg.jl index 12091091..db56f20a 100644 --- a/test/test_multilayerqg.jl +++ b/test/test_multilayerqg.jl @@ -172,7 +172,7 @@ function test_mqg_nonlinearadvection(dt, stepper, dev::Device=CPU(); nx, ny = 64, 66 Lx, Ly = 2π, 2π - gr = TwoDGrid(dev, nx, Lx, ny, Ly) + gr = TwoDGrid(dev; nx, Lx, ny, Ly) x, y = gridpoints(gr) k₀, l₀ = 2π/gr.Lx, 2π/gr.Ly # fundamental wavenumbers @@ -263,7 +263,7 @@ function test_mqg_linearadvection(dt, stepper, dev::Device=CPU(); nx, ny = 64, 66 Lx, Ly = 2π, 2π - gr = TwoDGrid(dev, nx, Lx, ny, Ly) + gr = TwoDGrid(dev; nx, Lx, ny, Ly) x, y = gridpoints(gr) k₀, l₀ = 2π/gr.Lx, 2π/gr.Ly # fundamental wavenumbers @@ -343,7 +343,7 @@ function test_mqg_energies(dev::Device=CPU(); dt=0.001, stepper="ForwardEuler", n=128, L=2π, nlayers=2, μ=0.0, ν=0.0, nν=1) nx, ny = 64, 66 Lx, Ly = 2π, 2π - gr = TwoDGrid(dev, nx, Lx, ny, Ly) + gr = TwoDGrid(dev; nx, Lx, ny, Ly) x, y = gridpoints(gr) k₀, l₀ = 2π/gr.Lx, 2π/gr.Ly # fundamental wavenumbers @@ -379,7 +379,7 @@ function test_mqg_energysinglelayer(dev::Device=CPU(); dt=0.001, stepper="ForwardEuler", nlayers=1, μ=0.0, ν=0.0, nν=1) nx, Lx = 64, 2π ny, Ly = 64, 3π - gr = TwoDGrid(dev, nx, Lx, ny, Ly) + gr = TwoDGrid(dev; nx, Lx, ny, Ly) x, y = gridpoints(gr) k₀, l₀ = 2π/gr.Lx, 2π/gr.Ly # fundamental wavenumbers @@ -407,7 +407,7 @@ initializing it with a flow field whose fluxes are known. function test_mqg_fluxes(dev::Device=CPU(); dt=0.001, stepper="ForwardEuler", n=128, L=2π, nlayers=2, μ=0.0, ν=0.0, nν=1) nx, ny = 128, 126 Lx, Ly = 2π, 2π - gr = TwoDGrid(dev, nx, Lx, ny, Ly) + gr = TwoDGrid(dev; nx, Lx, ny, Ly) x, y = gridpoints(gr) k₀, l₀ = 2π/gr.Lx, 2π/gr.Ly # fundamental wavenumbers @@ -447,7 +447,7 @@ function test_mqg_fluxessinglelayer(dev::Device=CPU(); dt=0.001, stepper="Forwar nx, ny = 128, 126 Lx, Ly = 2π, 2π - gr = TwoDGrid(dev, nx, Lx, ny, Ly) + gr = TwoDGrid(dev; nx, Lx, ny, Ly) x, y = gridpoints(gr) k₀, l₀ = 2π/gr.Lx, 2π/gr.Ly # fundamental wavenumbers @@ -475,7 +475,7 @@ given `q` or `ψ` respectively. function test_mqg_setqsetψ(dev::Device=CPU(); dt=0.001, stepper="ForwardEuler", n=64, L=2π, nlayers=2, μ=0.0, ν=0.0, nν=1) nx, ny = 32, 34 L = 2π - gr = TwoDGrid(dev, nx, L, ny, L) + gr = TwoDGrid(dev; nx, Lx=L, ny, Ly=L) x, y = gridpoints(gr) k₀, l₀ = 2π/gr.Lx, 2π/gr.Ly # fundamental wavenumbers @@ -521,7 +521,7 @@ Tests that `Params` constructor works with both mean flow `U` being a floats function test_mqg_paramsconstructor(dev::Device=CPU(); dt=0.001, stepper="ForwardEuler") nx, ny = 32, 34 L = 2π - gr = TwoDGrid(dev, nx, L, ny, L) + gr = TwoDGrid(dev; nx, Lx=L, ny, Ly=L) nlayers = 2 # these choice of parameters give the f₀, g = 1, 1 # desired PV-streamfunction relations diff --git a/test/test_singlelayerqg.jl b/test/test_singlelayerqg.jl index cba25c62..2cca28e3 100644 --- a/test/test_singlelayerqg.jl +++ b/test/test_singlelayerqg.jl @@ -309,7 +309,7 @@ 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) + grid = TwoDGrid(dev; nx, Lx, ny, Ly) k₀, l₀ = 2π/Lx, 2π/Ly # fundamental wavenumbers x, y = gridpoints(grid) @@ -340,7 +340,7 @@ Tests the kinetic and potential energy for an equivalent barotropic SingleLayerQ 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) + grid = TwoDGrid(dev; nx, Lx, ny, Ly) k₀, l₀ = 2π/Lx, 2π/Ly # fundamental wavenumbers x, y = gridpoints(grid) diff --git a/test/test_surfaceqg.jl b/test/test_surfaceqg.jl index 7de12de1..17d70e69 100644 --- a/test/test_surfaceqg.jl +++ b/test/test_surfaceqg.jl @@ -70,7 +70,7 @@ function test_sqg_kineticenergy_buoyancyvariance(dev::Device=CPU()) nx, Lx = 128, 2π ny, Ly = 128, 3π - grid = TwoDGrid(dev, nx, Lx, ny, Ly) + grid = TwoDGrid(dev; nx, Lx, ny, Ly) x, y = gridpoints(grid) k₀, l₀ = 2π/grid.Lx, 2π/grid.Ly # fundamental wavenumbers diff --git a/test/test_utils.jl b/test/test_utils.jl index 7724e1e0..fd566d91 100644 --- a/test/test_utils.jl +++ b/test/test_utils.jl @@ -11,13 +11,13 @@ function testpeakedisotropicspectrum(dev::Device=CPU()) ρtest = ρ[ (ρ.>15.0) .& (ρ.<=17.5)] qhρtest = qhρ[ (ρ.>15.0) .& (ρ.<=17.5)] - CUDA.@allowscalar isapprox(abs.(qhρtest)/abs(qhρtest[1]), (ρtest/ρtest[1]).^(-2), rtol=5e-3) + return CUDA.@allowscalar isapprox(abs.(qhρtest)/abs(qhρtest[1]), (ρtest/ρtest[1]).^(-2), rtol=5e-3) end -function testpeakedisotropicspectrum_rectangledomain() +function testpeakedisotropicspectrum_rectangledomain(dev::Device=CPU()) nx, ny = 32, 34 Lx, Ly = 2π, 3π - grid = TwoDGrid(nx, Lx, ny, Ly) + grid = TwoDGrid(dev; nx, Lx, ny, Ly) k0, E0 = 6, 0.5 qi = peakedisotropicspectrum(grid, k0, E0; allones=true) end From 7aadfdeb0a1e5f121e539a93de45540814fa0181 Mon Sep 17 00:00:00 2001 From: "Navid C. Constantinou" Date: Fri, 2 Sep 2022 14:04:51 +1000 Subject: [PATCH 6/7] update to latest grid refactor --- src/surfaceqg.jl | 2 +- src/twodnavierstokes.jl | 2 +- test/test_twodnavierstokes.jl | 6 +++--- 3 files changed, 5 insertions(+), 5 deletions(-) diff --git a/src/surfaceqg.jl b/src/surfaceqg.jl index e68e7439..722a3d9a 100644 --- a/src/surfaceqg.jl +++ b/src/surfaceqg.jl @@ -83,7 +83,7 @@ function Problem(dev::Device=CPU(); equation = Equation(params, grid) - return FourierFlows.Problem(equation, stepper, dt, grid, vars, params, dev) + return FourierFlows.Problem(equation, stepper, dt, grid, vars, params) end diff --git a/src/twodnavierstokes.jl b/src/twodnavierstokes.jl index bf9dbfd4..c3a2e918 100644 --- a/src/twodnavierstokes.jl +++ b/src/twodnavierstokes.jl @@ -91,7 +91,7 @@ function Problem(dev::Device=CPU(); equation = Equation(params, grid) - return FourierFlows.Problem(equation, stepper, dt, grid, vars, params, dev) + return FourierFlows.Problem(equation, stepper, dt, grid, vars, params) end diff --git a/test/test_twodnavierstokes.jl b/test/test_twodnavierstokes.jl index 209b2571..6c589aa1 100644 --- a/test/test_twodnavierstokes.jl +++ b/test/test_twodnavierstokes.jl @@ -246,7 +246,7 @@ function test_twodnavierstokes_advection(dt, stepper, dev::Device=CPU(); n=128, return nothing end - prob = TwoDNavierStokes.Problem(dev; nx=n, Lx=L, ν=ν, nν=nν, μ=μ, nμ=nμ, dt=dt, stepper=stepper, calcF=calcF!, stochastic=false) + prob = TwoDNavierStokes.Problem(dev; nx=n, Lx=L, ν, nν, μ, nμ, dt, stepper, calcF=calcF!, stochastic=false) TwoDNavierStokes.set_ζ!(prob, ζf) @@ -261,7 +261,7 @@ function test_twodnavierstokes_energyenstrophy(dev::Device=CPU()) nx, Lx = 128, 2π ny, Ly = 126, 3π - grid = TwoDGrid(dev, nx, Lx, ny, Ly) + grid = TwoDGrid(dev; nx, Lx, ny, Ly) x, y = gridpoints(grid) k₀, l₀ = 2π/grid.Lx, 2π/grid.Ly # fundamental wavenumbers @@ -271,7 +271,7 @@ function test_twodnavierstokes_energyenstrophy(dev::Device=CPU()) energy_calc = 29/9 enstrophy_calc = 2701/162 - prob = TwoDNavierStokes.Problem(dev; nx=nx, Lx=Lx, ny=ny, Ly=Ly, stepper="ForwardEuler") + prob = TwoDNavierStokes.Problem(dev; nx, Lx, ny, Ly, stepper="ForwardEuler") sol, cl, v, p, g = prob.sol, prob.clock, prob.vars, prob.params, prob.grid; From 4c8839e11ffa71b38ddf86b1dcc849e3bb32882e Mon Sep 17 00:00:00 2001 From: "Navid C. Constantinou" Date: Fri, 2 Sep 2022 14:13:11 +1000 Subject: [PATCH 7/7] bump minor release --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 31715bbd..3c607655 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 co-contributors"] -version = "0.14.2" +version = "0.15.0" [deps] CUDA = "052768ef-5323-5732-b1bb-66c8b64840ba"