Shear Flow with Periodic Boundaries
This simulation creates a periodic simulation with an ocean shear flow current. It is very similar to the example created in the tutorial. The tutorial explains the setup in greater detail. It is a good starter simulation for understanding how the code works and playing around with a basic setup.
This simulation shows how to setup a basic domain where all four boundaries are periodic. It also creates non-constant, shear ocean field. The ocean u-velocities are zero at the minimum and maximum y-extents (constant across x-values). The u values increase from the top and bottom of the domain towards the center to 0.5m/s, again constant across x-values, varying with y-values.
using Subzero, CairoMakie, GeoInterfaceMakie
using JLD2, Random, StatisticsUser Inputs
const FT = Float64 # Float type used to run simulation
const Lx = 1e5 # grid x-length
const Ly = 1e5 # grid y-length
const Δgrid = 2e3 # grid cell edge-size
const hmean = 0.25 # mean floe height
const Δh = 0.0 # difference in floe heights - here all floes are the same height
const Δt = 20 # timestep
const nΔt = 5000; # number of timesteps to runGrid Creation
grid = RegRectilinearGrid(; x0 = 0.0, xf = Lx, y0 = 0.0, yf = Ly, Δx = Δgrid, Δy = Δgrid)RegRectilinearGrid{Float64}
⊢x extent (0.0 to 100000.0) with 50 grid cells of size 2000.0 m
∟y extent (0.0 to 100000.0) with 50 grid cells of size 2000.0 mDomain Creation
nboundary = PeriodicBoundary(North; grid)
sboundary = PeriodicBoundary(South; grid)
eboundary = PeriodicBoundary(East; grid)
wboundary = PeriodicBoundary(West; grid)
domain = Domain(; north = nboundary, south = sboundary, east = eboundary, west = wboundary)Domain
⊢Northern boundary of type PeriodicBoundary{North, Float64}
⊢Southern boundary of type PeriodicBoundary{South, Float64}
⊢Eastern boundary of type PeriodicBoundary{East, Float64}
⊢Western boundary of type PeriodicBoundary{West, Float64}
∟0-element TopograpahyElement{Float64} listOcean Creation
half_nx = fld(grid.Nx, 2)
u_vec = [range(0, 0.5, length = half_nx); 0.5; range(0.5, 0, length = half_nx)]
u_row = transpose(u_vec)
uvels = repeat(u_row, outer = (grid.Ny + 1, 1))
ocean = Ocean(; u = uvels, v = 0, temp = 0, grid)Ocean{Float64}
⊢Vector fields of dimension (51, 51)
⊢Tracer fields of dimension (51, 51)
⊢Average u-velocity of: 0.2549 m/s
⊢Average v-velocity of: 0.0 m/s
∟Average temperature of: 0.0 CAtmosphere Creation
atmos = Atmos(; u = 0.0, v = 0.0, temp = -1.0, grid)Atmos{Float64}
⊢Vector fields of dimension (51, 51)
⊢Tracer fields of dimension (51, 51)
⊢Average u-velocity of: 0.0 m/s
⊢Average v-velocity of: 0.0 m/s
∟Average temperature of: -1.0 CFloe Creation
floe_settings = FloeSettings(subfloe_point_generator = SubGridPointsGenerator(; grid, npoint_per_cell = 2))
floe_generator = VoronoiTesselationFieldGenerator(; nfloes = 50, concentrations = [0.75], hmean, Δh)
floe_arr = initialize_floe_field(FT; generator = floe_generator, domain, rng = Xoshiro(1), floe_settings)51-element StructArray(::Vector{GeoInterface.Wrappers.Polygon{false, false, Vector{GeoInterface.Wrappers.LinearRing{false, false, Vector{Tuple{Float64, Float64}}, Nothing, Nothing}}, Nothing, Nothing}}, ::Vector{Vector{Float64}}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Vector{Float64}}, ::Vector{Vector{Float64}}, ::Vector{Vector{Float64}}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Subzero.Status}, ::Vector{Int64}, ::Vector{Int64}, ::Vector{Vector{Int64}}, ::Vector{Vector{Int64}}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Matrix{Float64}}, ::Vector{Float64}, ::Vector{Matrix{Float64}}, ::Vector{Int64}, ::Vector{Matrix{Float64}}, ::Vector{Matrix{Float64}}, ::Vector{Matrix{Float64}}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Float64}) with eltype Floe{Float64}:
Floe{Float64}
⊢Centroid of (34416.62809, 85023.74365) m
⊢Height of 0.25 m
⊢Area of 1.73342442274e7 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (82415.80025, 17595.45063) m
⊢Height of 0.25 m
⊢Area of 1.1417469701719e8 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (11656.4435, 93622.62301) m
⊢Height of 0.25 m
⊢Area of 2.1638856747361e8 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (31436.45811, 7814.63641) m
⊢Height of 0.25 m
⊢Area of 1.7917676347892e8 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (84666.43212, 93880.25523) m
⊢Height of 0.25 m
⊢Area of 1.696485653242e8 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (85872.18615, 50222.25842) m
⊢Height of 0.25 m
⊢Area of 2.0843316884575e8 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (71588.52548, 6806.77717) m
⊢Height of 0.25 m
⊢Area of 2.9573855564447e8 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (66409.08867, 50729.61018) m
⊢Height of 0.25 m
⊢Area of 1.4576742987025e8 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (81623.9572, 85028.5507) m
⊢Height of 0.25 m
⊢Area of 1.0964450041659e8 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (79809.59027, 40354.0117) m
⊢Height of 0.25 m
⊢Area of 2.117383616556e8 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
⋮
Floe{Float64}
⊢Centroid of (16488.30754, 65027.98406) m
⊢Height of 0.25 m
⊢Area of 1.9749028464881e8 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (32195.79981, 36824.15059) m
⊢Height of 0.25 m
⊢Area of 1.9236323809466e8 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (95773.67734, 46711.12264) m
⊢Height of 0.25 m
⊢Area of 1.7966044078233e8 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (33107.15163, 66820.16937) m
⊢Height of 0.25 m
⊢Area of 5.161377030363e7 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (11867.15435, 29800.9169) m
⊢Height of 0.25 m
⊢Area of 1.1646786042179e8 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (85002.91627, 76201.74791) m
⊢Height of 0.25 m
⊢Area of 6.679433210086e7 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (40771.33883, 54661.70948) m
⊢Height of 0.25 m
⊢Area of 2.4275388020636e8 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (75571.0654, 63500.51622) m
⊢Height of 0.25 m
⊢Area of 3.2777213000065e8 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Floe{Float64}
⊢Centroid of (5472.34891, 73418.50069) m
⊢Height of 0.25 m
⊢Area of 1.9384289442661e8 m^2
∟Velocity of (u, v, ξ) of (0.0, 0.0, 0.0) in (m/s, m/s, rad/s)
Model Creation
model = Model(; grid, ocean, atmos, domain, floes = floe_arr)Model{Float64, ...}
⊢RegRectilinearGrid{Float64}
⊢x extent (0.0 to 100000.0) with 50 grid cells of size 2000.0 m
∟y extent (0.0 to 100000.0) with 50 grid cells of size 2000.0 m
⊢Domain
⊢Northern boundary of type PeriodicBoundary{North, Float64}
⊢Southern boundary of type PeriodicBoundary{South, Float64}
⊢Eastern boundary of type PeriodicBoundary{East, Float64}
⊢Western boundary of type PeriodicBoundary{West, Float64}
∟0-element TopograpahyElement{Float64} list
⊢Ocean{Float64}
⊢Vector fields of dimension (51, 51)
⊢Tracer fields of dimension (51, 51)
⊢Average u-velocity of: 0.2549 m/s
⊢Average v-velocity of: 0.0 m/s
∟Average temperature of: 0.0 C
⊢Atmos{Float64}
⊢Vector fields of dimension (51, 51)
⊢Tracer fields of dimension (51, 51)
⊢Average u-velocity of: 0.0 m/s
⊢Average v-velocity of: 0.0 m/s
∟Average temperature of: -1.0 C
∟Floe List:
⊢Number of floes: 51
⊢Total floe area: 7.68058281616688e9
∟Average floe height: 0.25Constants Creation
modulus = 1.5e3*(mean(sqrt.(floe_arr.area)) + minimum(sqrt.(floe_arr.area)))
consts = Constants(; E = modulus)Constants{Float64}
⊢ ρo = 1027.0
⊢ ρa = 1.2
⊢ Cd_io = 0.003
⊢ Cd_ia = 0.001
⊢ Cd_ao = 0.00125
⊢ f = 0.00014
⊢ turnθ = 0.2617993877991494
⊢ L = 293000.0
⊢ k = 2.14
⊢ ν = 0.3
⊢ μ = 0.2
⊢ E = 2.065436905321266e7
Output Creation
dir = "shear_flow"
init_fn, floe_fn = "shear_flow_init_state.jld2", "shear_flow_floes.jld2"
initwriter = InitialStateOutputWriter(dir = dir, filename = init_fn, overwrite = true)
floewriter = FloeOutputWriter(50, dir = dir, filename = floe_fn, overwrite = true)
writers = OutputWriters(initwriter, floewriter)OutputWriters
⊢1 InitialStateOuputWriter(s)
⊢0 CheckpointOutputWriter(s)
⊢1 FloeOutputWriter(s)
∟0 GridOutputWriter(s)Simulation Creation
simulation = Simulation(; model, consts, writers, Δt, nΔt, floe_settings, verbose = true, rng = Xoshiro(1))Simulation
⊢Timestep: 20 seconds
⊢Runtime: 5000 timesteps
⊢RNG: Xoshiro(0xfff0241072ddab67, 0xc53bc12f4c3f0b4e, 0x56d451780b2dd4ba, 0x50a4aa153d208dd8, 0x3649a58b3b63d5db)
⊢verbose: true
⊢model
⊢consts
⊢floe_settings
⊢collision_settings
∟ ...Running the Simulation
@time run!(simulation)sim is running!
0 timesteps
50 timesteps
100 timesteps
150 timesteps
200 timesteps
250 timesteps
300 timesteps
350 timesteps
400 timesteps
450 timesteps
500 timesteps
550 timesteps
600 timesteps
650 timesteps
700 timesteps
750 timesteps
800 timesteps
850 timesteps
900 timesteps
950 timesteps
1000 timesteps
1050 timesteps
1100 timesteps
1150 timesteps
1200 timesteps
1250 timesteps
1300 timesteps
1350 timesteps
1400 timesteps
1450 timesteps
1500 timesteps
1550 timesteps
1600 timesteps
1650 timesteps
1700 timesteps
1750 timesteps
1800 timesteps
1850 timesteps
1900 timesteps
1950 timesteps
2000 timesteps
2050 timesteps
2100 timesteps
2150 timesteps
2200 timesteps
2250 timesteps
2300 timesteps
2350 timesteps
2400 timesteps
2450 timesteps
2500 timesteps
2550 timesteps
2600 timesteps
2650 timesteps
2700 timesteps
2750 timesteps
2800 timesteps
2850 timesteps
2900 timesteps
2950 timesteps
3000 timesteps
3050 timesteps
3100 timesteps
3150 timesteps
3200 timesteps
3250 timesteps
3300 timesteps
3350 timesteps
3400 timesteps
3450 timesteps
3500 timesteps
3550 timesteps
3600 timesteps
3650 timesteps
3700 timesteps
3750 timesteps
3800 timesteps
3850 timesteps
3900 timesteps
3950 timesteps
4000 timesteps
4050 timesteps
4100 timesteps
4150 timesteps
4200 timesteps
4250 timesteps
4300 timesteps
4350 timesteps
4400 timesteps
4450 timesteps
4500 timesteps
4550 timesteps
4600 timesteps
4650 timesteps
4700 timesteps
4750 timesteps
4800 timesteps
4850 timesteps
4900 timesteps
4950 timesteps
5000 timesteps
sim done running!
11.587965 seconds (71.44 M allocations: 6.322 GiB, 8.71% gc time, 0.29% compilation time)Plotting the Simulation
plot_sim(joinpath(dir, floe_fn), joinpath(dir, init_fn), Δt, joinpath(dir, "shear_flow.mp4"))Note that this is just using the built-in basic plotting. However, it is easy to write your own plotting code. See the source code for a basic outline.
This page was generated using Literate.jl.