Initial value problem

[maeckha]

1. Set up an initial value problem without particles, including visualization. Inspect the flow, tweak parameters, and develop an intuition for how changes to parameters change the resulting flow. Settle on a setup that you think will produce something sensible.
2. Add a very simple implementation of Lagrangian particles that are passively advected, with no “dynamics” / rock papers scissor rule. Set up visualization that combines the flow visuzliation with particle trajectors. Iterate on the setup so that the visualization is helpful and informative.
3. Implement incremental changes to the flow or particle dynamics that get you closer to the problem you’d like to simulate ultimately.

Per the approach above I have started with a uniform condition.

grid = RegularRectilinearGrid(size = (16, 16, 1), extent = (2π, 2π, 1))
model = NonhydrostaticModel(grid = grid)
set!(model, u=0.1, v=1.5)

What would the initial condition have to be, if a laminar flow was to enter left and leave right?

[glwagner]

Try

grid = RegularRectilinearGrid(size = (128, 128), x=(-5, 5), y=(-5, 5), topology=(Periodic, Periodic, Flat))
model = NonhydrostaticModel(grid = grid, tracers=:c, buoyancy=nothing)

initial_c(x, y, z) = exp(-x^2 - y^2)
set!(model, u=1, c=initial_c)

this initializes a uniform velocity from left to right, with a Gaussian tracer profile:

using UnicodePlots
heatmap(interior(model.tracers.c)[:, :, 1]')

producing

stepping the model forward shows the effect of advection by a uniform flow on the tracer distribution:

simulation = Simulation(model, Δt=1e-2, stop_iteration=100)
run!(simulation)
heatmap(interior(model.tracers.c)[:, :, 1]')

[glwagner]

This part of the code:

simulation = Simulation(model, Δt=1e-2, stop_iteration=100)
run!(simulation)

causes the model to step forward 100 time-steps using the QuasiAdamsBashforth2TimeStepper. As a result of this time-stepping, the tracer distribution moves from left to right. Try changing these two lines to

simulation = Simulation(model, Δt=1e-2, stop_iteration=400)
run!(simulation)

(increasing stop_iteration from 100 to 400) and see what happens.

[glwagner]

I believe after the particles are created, you have to add them to the model by writing

model = NonhydrostaticModel(grid = grid, tracers=:c, buoyancy=nothing, lagrangian_particles=lagrangian_particles)

I would also use

z₀ = zeros(n_particles)

[maeckha]

I have set the timesteps too small and therefore could not see enough progress in the heatmap. Stepping from 100 to 400 looks better. Don't know why I have stuck to incremental stepping for that long, thanks.

Does the 'u' here mean velocity?
set!(model, u=1, c=initial_c)

[glwagner]

On a rectilinear grid, u, v, w are the velocity components in the x, y, z direction.

[maeckha]

I believe after the particles are created, you have to add them to the model by writing

model = NonhydrostaticModel(grid = grid, tracers=:c, buoyancy=nothing, lagrangian_particles=lagrangian_particles)

I would also use

z₀ = zeros(n_particles)

model = NonhydrostaticModel(grid = grid, tracers=:c, buoyancy=nothing, particles=lagrangian_particles)