-
Notifications
You must be signed in to change notification settings - Fork 1
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
* Redo README text and make new files for other documentation sections * Update readme to make it nicer * Create title simulation * Change title video to gif * Add title gif * Remove spacing in readme
- Loading branch information
Showing
7 changed files
with
214 additions
and
135 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,97 @@ | ||
using Subzero, JLD2, CairoMakie, GeoInterfaceMakie, Random, Statistics | ||
|
||
const FT = Float64 | ||
const Lx = 8e4 | ||
const Ly = 8e4 | ||
const Δgrid = 2e3 | ||
const hmean = 0.25 | ||
const Δh = 0.0 | ||
const Δt = 20 | ||
const nΔt = 5000 | ||
|
||
# Model instantiation | ||
grid = RegRectilinearGrid((0.0, Lx), (0.0, Ly), Δgrid, Δgrid) | ||
|
||
uvels = repeat( | ||
[range(0, 0.2, length = 21); range(0.2, 0, length = 20)], | ||
outer = (1, 41), | ||
) | ||
ocean = Ocean( | ||
uvels', | ||
zeros(grid.Nx + 1, grid.Ny + 1), | ||
zeros(grid.Nx + 1, grid.Ny + 1), | ||
) | ||
atmos = Atmos(grid, 0.0, 0.0, -1.0) | ||
|
||
# Domain creation | ||
nboundary = PeriodicBoundary(North, grid) | ||
sboundary = PeriodicBoundary(South, grid) | ||
eboundary = PeriodicBoundary(East, grid) | ||
wboundary = PeriodicBoundary(West, grid) | ||
domain = Domain(nboundary, sboundary, eboundary, wboundary) | ||
|
||
# Floe creation | ||
floe_arr = initialize_floe_field(FT, 400, [0.85], domain, hmean, Δh; rng = Xoshiro(1)) | ||
|
||
# Model creation | ||
model = Model(grid, ocean, atmos, domain, floe_arr) | ||
|
||
# Simulation setup | ||
modulus = 1.5e3*(mean(sqrt.(floe_arr.area)) + minimum(sqrt.(floe_arr.area))) | ||
consts = Constants(E = modulus) | ||
|
||
# Run simulation | ||
dir = "output/title/" | ||
|
||
# Output setup | ||
floewriter = FloeOutputWriter(50, dir = dir, overwrite = true) | ||
writers = OutputWriters(floewriter) | ||
|
||
simulation = Simulation(; model, consts, Δt, nΔt, writers, verbose = true, rng = Xoshiro(1)) | ||
run!(simulation) | ||
|
||
function plot_logo(floe_fn, Lx, Ly, dir) | ||
# Floe Information | ||
file = jldopen(floe_fn) | ||
sim_polys = file["poly"] | ||
timesteps = keys(sim_polys) | ||
|
||
# Set up observables for recording (updated whenever `time[]` is set to a new value) | ||
time = Observable(timesteps[1]) # note these are strings as we index into a JLD2 file | ||
time_polys = @lift(sim_polys[$time]) | ||
|
||
# Set up figure | ||
fig = Figure(; size = (1200, 300)) | ||
|
||
floe_rgb = (165, 222, 242) | ||
floe_color = RGBf((floe_rgb ./ 255)...) | ||
ocean_rgb = (4, 31, 74) | ||
ocean_color = RGBf((ocean_rgb ./ 255)...) | ||
text_outline_rgb = (219, 76, 0) | ||
text_outline_color = RGBf((text_outline_rgb ./ 255)...) | ||
|
||
ax = Axis(fig[1, 1]; limits = (0.0, Lx, 0.0, Ly), backgroundcolor = ocean_color) | ||
hidedecorations!(ax) | ||
|
||
# Plot starting state (floes + topography) | ||
poly!(time_polys, color = floe_color, strokecolor = :black, strokewidth = 0.5) # floes | ||
text!( | ||
Lx/2, Ly/2; | ||
align = (:center, :center), | ||
text = "Subzero.jl", | ||
font = :bold, | ||
fontsize = 175, | ||
strokecolor = text_outline_color, | ||
strokewidth = 5, | ||
) | ||
|
||
# Create movie | ||
record(fig, dir*"title.gif", timesteps; framerate = 20) do t | ||
time[] = t # updating the time auto updates the related time_polys | ||
end | ||
|
||
# End movie and close file | ||
close(file) | ||
end | ||
|
||
plot_logo(dir*"floes.jld2", Lx, Ly, dir) |
Empty file.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -1,3 +1,5 @@ | ||
# Introduction | ||
|
||
Sea-ice model to explicitly simulate individual floe life cycles using complex discrete elements with time-evolving shapes. Easy to couple with Oceananigans for explorations of coupled ice-ocean dynamics. Ported from MATLAB to Julia. | ||
|
||
## Main concepts |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,78 @@ | ||
# Let’s run a basic simulation with initially stationary floes pushed into a collision boundary by a uniform, zonally flowing ocean. In this simulation, collisions between floes are on by default and we will enable floe fracturing. | ||
|
||
# Create environment | ||
grid = RegRectilinearGrid( | ||
(0, 1e5), # xbounds | ||
(0, 1e5), # ybounds | ||
1e4, # grid cell width | ||
1e4, # grid cell height | ||
) | ||
ocean = Ocean(grid, 0.25, 0.0, 0.0) # 0.25m/s u-velocity, 0m/s v-velocity, 0C temperature in all grid cells | ||
atmos = Atmos(grid, 0.0, 0.1, -1.0) # 0m/s u-velocity, 0.1m/s v-velocity, -1C temperature in all grid cells | ||
# Create domain | ||
domain = Domain( | ||
CollisionBoundary(North, grid), | ||
CollisionBoundary(South, grid), | ||
CollisionBoundary(East, grid), | ||
CollisionBoundary(West, grid), | ||
) | ||
# Create floes | ||
floe_settings = FloeSettings(stress_calculator = DecayAreaScaledCalculator(), min_floe_area = 1e5) | ||
floe_arr = initialize_floe_field( | ||
Float64, | ||
100, # number of floes | ||
[0.7], # floe concentration | ||
domain, | ||
0.5, # average floe height | ||
0.05; # floe height variability | ||
floe_settings = floe_settings, | ||
) | ||
# Create model | ||
model = Model(grid, ocean, atmos, domain, floe_arr) | ||
|
||
# Create settings | ||
modulus = 1.5e3*(mean(sqrt.(floe_arr.area)) + minimum(sqrt.(floe_arr.area))) | ||
consts = Constants(E = modulus) | ||
|
||
fracture_settings = FractureSettings( | ||
fractures_on = true, | ||
criteria = HiblerYieldCurve(floe_arr), | ||
Δt = 75, | ||
npieces = 3, | ||
deform_on = true, | ||
) | ||
# Create output writers | ||
dir = "output/sim" | ||
initwriter = InitialStateOutputWriter( | ||
dir = dir, | ||
filename = "initial_state.jld2", | ||
overwrite = true, | ||
) | ||
floewriter = FloeOutputWriter( | ||
100, | ||
dir = dir, | ||
filename = "floes.jld2", | ||
overwrite = true, | ||
) | ||
writers = OutputWriters(initwriter, floewriter) | ||
# Create simulation | ||
Δt = 10 | ||
simulation = Simulation( | ||
model = model, | ||
consts = consts, | ||
Δt = Δt, # 10 second timesteps | ||
nΔt = 10000, # Run for 10000 timesteps | ||
verbose = true; | ||
fracture_settings = fracture_settings, | ||
writers = writers, | ||
) | ||
# Run simulation | ||
run!(simulation) | ||
|
||
# Plot simulation | ||
plot_sim( | ||
"output/sim/floes.jld2", | ||
"output/sim/initial_state.jld2", | ||
Δt, | ||
"output/sim/example.mp4", | ||
) |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters