Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

[WIP] Simplify bench tests #10

Open
wants to merge 3 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
64 changes: 29 additions & 35 deletions GeophysicalFlows.jl/GeophysicalFlows_decaying2Dturbulence.jl
Original file line number Diff line number Diff line change
@@ -1,70 +1,64 @@
using FourierFlows, Printf, Random, FFTW

using Random: seed!
using FourierFlows, Printf

import GeophysicalFlows.TwoDNavierStokes
import GeophysicalFlows.TwoDNavierStokes: energy, enstrophy
import GeophysicalFlows: peakedisotropicspectrum

using Random: seed!
using Statistics: mean

dev = CPU() # Device (CPU/GPU)

# Parameters
n = 256
L = 2π
= 2
ν = 0.0
dt = 1e-3
nstepsjit = 10 # call stepforward for nstepsjit to force compilation
n = 256
L = 2π
ν = 1e-4
nν = 1
dt = 1e-3
nsteps_jit = 10 # call stepforward for nsteps_jit to force compilation
nsteps = 5000
nothingfunction(args...) = nothing


gr = TwoDGrid(dev, n, L, n, L; nthreads=1)
pr = TwoDNavierStokes.Params(ν, nν)
vs = TwoDNavierStokes.Vars(dev, gr)
eq = TwoDNavierStokes.Equation(pr, gr)

prob = FourierFlows.Problem(eq, "FilteredAB3", dt, gr, vs, pr, dev)
filter = FourierFlows.makefilter(prob.eqn)
stepper = "AB3"
floattype = Float64

# some aliases
sol, cl, vs, gr = prob.sol, prob.clock, prob.vars, prob.grid
grid2D = TwoDGrid(dev, n, L, n, L; nthreads=1, T=floattype)
params = TwoDNavierStokes.Params(ν, nν)
vars = TwoDNavierStokes.Vars(dev, grid2D)
equation = TwoDNavierStokes.Equation(params, grid2D)

x, y = gridpoints(gr)
prob = FourierFlows.Problem(equation, stepper, dt, grid2D, vars, params, dev)

# Initial condition closely following pyqg barotropic example
# that reproduces the results of the paper by McWilliams (1984)
# Random initial condition
seed!(1234)
k0, E0 = 6, 0.5
zetai = peakedisotropicspectrum(gr, k0, E0, mask=filter)
TwoDNavierStokes.set_zeta!(prob, zetai)
zeta_initial = randn(floattype, (n, n))
zeta_initial = zeta_initial .- mean(zeta_initial) # make sure initial condition has zero mean
TwoDNavierStokes.set_zeta!(prob, zeta_initial)

for j=1:nstepsjit #just-in-time compilation
for j=1:nsteps_jit #just-in-time compilation
stepforward!(prob)
end

TwoDNavierStokes.set_zeta!(prob, zetai)
TwoDNavierStokes.set_zeta!(prob, zeta_initial)

startwalltime = time()
while cl.step < nsteps
while prob.clock.step < nsteps
stepforward!(prob)
dealias!(prob.sol, prob.grid)
end

println(round((time()-startwalltime)/(nsteps-1)*1000, digits=3), " ms per time-step")

# using PyPlot
# x, y = gridpoints(prob.grid)
# TwoDNavierStokes.updatevars!(prob)
# figure(figsize=(10, 4))
# subplot(121)
# pcolormesh(x, y, zetai)
# pcolormesh(x, y, zeta_initial)
# xlabel("x")
# ylabel("y")
# title("vorticity @ t=0")
# axis("square")
# subplot(122)
# pcolormesh(x, y, vs.zeta)
# pcolormesh(x, y, prob.vars.zeta)
# axis("square")
# xlabel("x")
# ylabel("y")
# title("vorticity @ t="*string(round(cl.t, digits=2)))
# title("vorticity @ t="*string(round(prob.clock.t, digits=2)))
# savefig("GeophysicalFlows_n"*string(n)*".png", dpi=400)
Binary file modified GeophysicalFlows.jl/GeophysicalFlows_n256.png
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Loading