Sea ice advected by an atmospheric anticyclone

using Oceananigans
using Oceananigans.Units
using Oceananigans.Operators
using ClimaSeaIce
using Printf
using ClimaSeaIce.SeaIceMomentumEquations
using ClimaSeaIce.Rheologies

The experiment found in the paper: Simulating Linear Kinematic Features in Viscous-Plastic Sea Ice Models on Quadrilateral and Triangular Grids With Different Variable Staggering

arch = CPU()

L  = 512kilometers
𝓋ₒ = 0.01 # m / s maximum ocean speed
𝓋ₐ = 30.0 # m / s maximum atmospheric speed modifier
Cᴰ = 1.2e-3 # Atmosphere - sea ice drag coefficient
ρₐ = 1.3  # kg/m³

1.3

2 km domain

grid = RectilinearGrid(arch;
                       size = (128, 128),
                          x = (0, L),
                          y = (0, L),
                   topology = (Bounded, Bounded, Flat))

#####
##### Value boundary conditions for velocities
#####

u_bcs = FieldBoundaryConditions(north = ValueBoundaryCondition(0),
                                south = ValueBoundaryCondition(0))

v_bcs = FieldBoundaryConditions(west = ValueBoundaryCondition(0),
                                east = ValueBoundaryCondition(0))

#####
##### Ocean sea-ice stress
#####

Oceananigans.FieldBoundaryConditions, with boundary conditions
├── west: ValueBoundaryCondition: 0
├── east: ValueBoundaryCondition: 0
├── south: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── north: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── bottom: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── top: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
└── immersed: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)

Constant ocean velocities corresponding to a cyclonic eddy

Uₒ = XFaceField(grid)
Vₒ = YFaceField(grid)

set!(Uₒ, (x, y) -> 𝓋ₒ * (2y - L) / L)
set!(Vₒ, (x, y) -> 𝓋ₒ * (L - 2x) / L)

Oceananigans.BoundaryConditions.fill_halo_regions!(Uₒ)
Oceananigans.BoundaryConditions.fill_halo_regions!(Vₒ)

struct ExplicitOceanSeaIceStress{U, V, C}
    u    :: U
    v    :: V
    ρₒCᴰ :: C
end

We extend the τx and τy methods to compute the time-dependent stress

import ClimaSeaIce.SeaIceMomentumEquations: τx, τy

@inline function τx(i, j, k, grid, τ::ExplicitOceanSeaIceStress, clock, fields)
    Δu = @inbounds fields.u[i, j, k] - τ.u[i, j, k]
    Δv = ℑxyᶠᶜᵃ(i, j, k, grid, τ.v) - ℑxyᶠᶜᵃ(i, j, k, grid, fields.v)
    return - τ.ρₒCᴰ * sqrt(Δu^2 + Δv^2) * Δu
end

@inline function τy(i, j, k, grid, τ::ExplicitOceanSeaIceStress, clock, fields)
    Δu = ℑxyᶜᶠᵃ(i, j, k, grid, τ.u) - ℑxyᶜᶠᵃ(i, j, k, grid, fields.u)
    Δv = @inbounds fields.v[i, j, k] - τ.v[i, j, k]
    return - τ.ρₒCᴰ * sqrt(Δu^2 + Δv^2) * Δv
end

τᵤₒ = τᵥₒ = ExplicitOceanSeaIceStress(Uₒ, Vₒ, 5.5)

####
#### Atmosphere - sea ice stress
####

Uₐ = XFaceField(grid)
Vₐ = YFaceField(grid)

128×129×1 Field{Oceananigans.Grids.Center, Oceananigans.Grids.Face, Oceananigans.Grids.Center} on Oceananigans.Grids.RectilinearGrid on Oceananigans.Architectures.CPU
├── grid: 128×128×1 RectilinearGrid{Float64, Oceananigans.Grids.Bounded, Oceananigans.Grids.Bounded, Oceananigans.Grids.Flat} on Oceananigans.Architectures.CPU with 3×3×0 halo
├── boundary conditions: FieldBoundaryConditions
│   └── west: ZeroFlux, east: ZeroFlux, south: Nothing, north: Nothing, bottom: Nothing, top: Nothing, immersed: ZeroFlux
└── data: 134×135×1 OffsetArray(::Array{Float64, 3}, -2:131, -2:132, 1:1) with eltype Float64 with indices -2:131×-2:132×1:1
    └── max=0.0, min=0.0, mean=0.0

Atmosphere - sea ice stress

τᵤₐ = Field(ρₐ * Cᴰ * sqrt(Uₐ^2 + Vₐ^2) * Uₐ)
τᵥₐ = Field(ρₐ * Cᴰ * sqrt(Uₐ^2 + Vₐ^2) * Vₐ)

129×128×1 Field{Oceananigans.Grids.Face, Oceananigans.Grids.Center, Oceananigans.Grids.Center} on Oceananigans.Grids.RectilinearGrid on Oceananigans.Architectures.CPU
├── grid: 128×128×1 RectilinearGrid{Float64, Oceananigans.Grids.Bounded, Oceananigans.Grids.Bounded, Oceananigans.Grids.Flat} on Oceananigans.Architectures.CPU with 3×3×0 halo
├── boundary conditions: FieldBoundaryConditions
│   └── west: Nothing, east: Nothing, south: ZeroFlux, north: ZeroFlux, bottom: Nothing, top: Nothing, immersed: ZeroFlux
├── operand: BinaryOperation at (Face, Center, Center)
├── status: time=0.0
└── data: 135×134×1 OffsetArray(::Array{Float64, 3}, -2:132, -2:131, 1:1) with eltype Float64 with indices -2:132×-2:131×1:1
    └── max=0.0, min=0.0, mean=0.0

Atmospheric velocities corresponding to an anticyclonic eddy moving north-east

@inline center(t) = 256kilometers + 51.2kilometers * t / 86400
@inline radius(x, y, t)  = sqrt((x - center(t))^2 + (y - center(t))^2)
@inline speed(x, y, t)   = 1 / 100 * exp(- radius(x, y, t) / 100kilometers)

@inline ua_time(x, y, t) = - 𝓋ₐ * speed(x, y, t) * (cosd(72) * (x - center(t)) + sind(72) * (y - center(t))) / 1000
@inline va_time(x, y, t) = - 𝓋ₐ * speed(x, y, t) * (cosd(72) * (y - center(t)) - sind(72) * (x - center(t))) / 1000

va_time (generic function with 1 method)

Initialize the stress at time t = 0

set!(Uₐ, (x, y) -> ua_time(x, y, 0))
set!(Vₐ, (x, y) -> va_time(x, y, 0))
compute!(τᵤₐ)
compute!(τᵥₐ)

#####
##### Numerical details
#####

129×128×1 Field{Oceananigans.Grids.Face, Oceananigans.Grids.Center, Oceananigans.Grids.Center} on Oceananigans.Grids.RectilinearGrid on Oceananigans.Architectures.CPU
├── grid: 128×128×1 RectilinearGrid{Float64, Oceananigans.Grids.Bounded, Oceananigans.Grids.Bounded, Oceananigans.Grids.Flat} on Oceananigans.Architectures.CPU with 3×3×0 halo
├── boundary conditions: FieldBoundaryConditions
│   └── west: Nothing, east: Nothing, south: ZeroFlux, north: ZeroFlux, bottom: Nothing, top: Nothing, immersed: ZeroFlux
├── operand: BinaryOperation at (Face, Center, Center)
├── status: time=0.0
└── data: 135×134×1 OffsetArray(::Array{Float64, 3}, -2:132, -2:131, 1:1) with eltype Float64 with indices -2:132×-2:131×1:1
    └── max=0.189842, min=-0.189842, mean=-6.21814e-19

We use an elasto-visco-plastic rheology and WENO seventh order for advection of h and ℵ

momentum_equations = SeaIceMomentumEquation(grid;
                                            coriolis = FPlane(f = 1e-4),
                                            rheology = ViscousRheology(ν = 1000.0))
advection = WENO(; order = 7)

u_bcs = FieldBoundaryConditions(north = ValueBoundaryCondition(0),
                                south = ValueBoundaryCondition(0))

v_bcs = FieldBoundaryConditions(west = ValueBoundaryCondition(0),
                                east = ValueBoundaryCondition(0))

Oceananigans.FieldBoundaryConditions, with boundary conditions
├── west: ValueBoundaryCondition: 0
├── east: ValueBoundaryCondition: 0
├── south: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── north: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── bottom: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── top: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
└── immersed: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)

Define the model!

model = SeaIceModel(grid;
                    top_momentum_stress = (u = τᵤₐ, v = τᵥₐ),
                    bottom_momentum_stress = (u = τᵤₒ, v = τᵥₒ),
                    ice_dynamics = momentum_equations,
                    ice_thermodynamics = nothing, # No thermodynamics here
                    advection,
                    boundary_conditions = (u = u_bcs, v = v_bcs))

SeaIceModel{Oceananigans.Architectures.CPU, Oceananigans.Grids.RectilinearGrid}(time = 0 seconds, iteration = 0)
├── grid: 128×128×1 RectilinearGrid{Float64, Oceananigans.Grids.Bounded, Oceananigans.Grids.Bounded, Oceananigans.Grids.Flat} on Oceananigans.Architectures.CPU with 3×3×0 halo
├── ice_thermodynamics: Nothing
├── advection: WENO(order=7)
└── external_heat_fluxes: 
    ├── top: Nothing
    └── bottom: 0

Initial height field with perturbations around 0.3 m

h₀(x, y) = 0.3 + 0.005 * (sin(60 * x / 1000kilometers) + sin(30 * y / 1000kilometers))

h₀ (generic function with 1 method)

We start with a concentration of ℵ = 1

set!(model, h = h₀)
set!(model, ℵ = 1)

#####
##### Setup the simulation
#####

run the model for 2 days

simulation = Simulation(model, Δt = 10seconds, stop_time = 2hours)

Simulation of SeaIceModel
├── Next time step: 10 seconds
├── Elapsed wall time: 0 seconds
├── Wall time per iteration: NaN days
├── Stop time: 2 hours
├── Stop iteration: Inf
├── Wall time limit: Inf
├── Minimum relative step: 0.0
├── Callbacks: OrderedDict with 3 entries:
│   ├── stop_time_exceeded => Callback of stop_time_exceeded on IterationInterval(1)
│   ├── stop_iteration_exceeded => Callback of stop_iteration_exceeded on IterationInterval(1)
│   └── wall_time_limit_exceeded => Callback of wall_time_limit_exceeded on IterationInterval(1)
├── Output writers: OrderedDict with no entries
└── Diagnostics: OrderedDict with no entries

Remember to evolve the wind stress field in time!

function compute_wind_stress(sim)
    time = sim.model.clock.time
    @inline ua(x, y) = ua_time(x, y, time)
    @inline va(x, y) = va_time(x, y, time)
    set!(Uₐ, ua)
    set!(Vₐ, va)

    compute!(τᵤₐ)
    compute!(τᵥₐ)

    return nothing
end

simulation.callbacks[:top_stress] = Callback(compute_wind_stress, IterationInterval(1))

Callback of compute_wind_stress on IterationInterval(1)

Container to hold the data

htimeseries = []
ℵtimeseries = []
utimeseries = []
vtimeseries = []

# Callback function to collect the data from the `sim`ulation
function accumulate_timeseries(sim)
    h = sim.model.ice_thickness
    ℵ = sim.model.ice_concentration
    u = sim.model.velocities.u
    v = sim.model.velocities.v
    push!(htimeseries, deepcopy(interior(h)))
    push!(ℵtimeseries, deepcopy(interior(ℵ)))
    push!(utimeseries, deepcopy(interior(u)))
    push!(vtimeseries, deepcopy(interior(v)))
end

wall_time = [time_ns()]

function progress(sim)
    h = sim.model.ice_thickness
    ℵ = sim.model.ice_concentration
    u = sim.model.velocities.u
    v = sim.model.velocities.v

    hmax = maximum(interior(h))
    ℵmin = minimum(interior(ℵ))
    umax = maximum(interior(u)), maximum(interior(v))
    step_time = 1e-9 * (time_ns() - wall_time[1])

    @info @sprintf("Time: %s, Iteration %d, Δt %s, max(vel): (%.2e, %.2e), max(h): %.2f, min(ℵ): %.2f, wtime: %s \n",
                   prettytime(sim.model.clock.time),
                   sim.model.clock.iteration,
                   prettytime(sim.Δt),
                   umax..., hmax, ℵmin, prettytime(step_time))

     wall_time[1] = time_ns()
end

simulation.callbacks[:progress] = Callback(progress, IterationInterval(5))
simulation.callbacks[:save]     = Callback(accumulate_timeseries, IterationInterval(5))

run!(simulation)

using CairoMakie

[ Info: Initializing simulation...
[ Info: Time: 0 seconds, Iteration 0, Δt 10 seconds, max(vel): (0.00e+00, 0.00e+00), max(h): 0.31, min(ℵ): 1.00, wtime: 40.230 seconds
[ Info:     ... simulation initialization complete (6.618 seconds)
[ Info: Executing initial time step...
[ Info:     ... initial time step complete (8.453 seconds).
[ Info: Time: 50.000 seconds, Iteration 15, Δt 10 seconds, max(vel): (3.48e-02, 3.40e-02), max(h): 0.31, min(ℵ): 1.00, wtime: 9.163 seconds
[ Info: Time: 1.667 minutes, Iteration 30, Δt 10 seconds, max(vel): (6.70e-02, 6.56e-02), max(h): 0.31, min(ℵ): 1.00, wtime: 747.439 ms
[ Info: Time: 2.500 minutes, Iteration 45, Δt 10 seconds, max(vel): (9.48e-02, 9.31e-02), max(h): 0.31, min(ℵ): 1.00, wtime: 756.945 ms
[ Info: Time: 3.333 minutes, Iteration 60, Δt 10 seconds, max(vel): (1.18e-01, 1.16e-01), max(h): 0.31, min(ℵ): 1.00, wtime: 760.157 ms
[ Info: Time: 4.167 minutes, Iteration 75, Δt 10 seconds, max(vel): (1.35e-01, 1.34e-01), max(h): 0.31, min(ℵ): 1.00, wtime: 760.989 ms
[ Info: Time: 5 minutes, Iteration 90, Δt 10 seconds, max(vel): (1.49e-01, 1.47e-01), max(h): 0.31, min(ℵ): 1.00, wtime: 745.186 ms
[ Info: Time: 5.833 minutes, Iteration 105, Δt 10 seconds, max(vel): (1.59e-01, 1.57e-01), max(h): 0.31, min(ℵ): 1.00, wtime: 769.775 ms
[ Info: Time: 6.667 minutes, Iteration 120, Δt 10 seconds, max(vel): (1.66e-01, 1.65e-01), max(h): 0.31, min(ℵ): 1.00, wtime: 765.114 ms
[ Info: Time: 7.500 minutes, Iteration 135, Δt 10 seconds, max(vel): (1.71e-01, 1.70e-01), max(h): 0.31, min(ℵ): 1.00, wtime: 744.417 ms
[ Info: Time: 8.333 minutes, Iteration 150, Δt 10 seconds, max(vel): (1.74e-01, 1.73e-01), max(h): 0.31, min(ℵ): 1.00, wtime: 738.986 ms
[ Info: Time: 9.167 minutes, Iteration 165, Δt 10 seconds, max(vel): (1.77e-01, 1.76e-01), max(h): 0.31, min(ℵ): 1.00, wtime: 769.873 ms
[ Info: Time: 10.000 minutes, Iteration 180, Δt 10 seconds, max(vel): (1.78e-01, 1.78e-01), max(h): 0.31, min(ℵ): 1.00, wtime: 772.543 ms
[ Info: Time: 10.833 minutes, Iteration 195, Δt 10 seconds, max(vel): (1.80e-01, 1.79e-01), max(h): 0.31, min(ℵ): 1.00, wtime: 776.562 ms
[ Info: Time: 11.667 minutes, Iteration 210, Δt 10 seconds, max(vel): (1.80e-01, 1.80e-01), max(h): 0.31, min(ℵ): 1.00, wtime: 778.905 ms
[ Info: Time: 12.500 minutes, Iteration 225, Δt 10 seconds, max(vel): (1.81e-01, 1.81e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 790.037 ms
[ Info: Time: 13.333 minutes, Iteration 240, Δt 10 seconds, max(vel): (1.81e-01, 1.81e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 797.758 ms
[ Info: Time: 14.167 minutes, Iteration 255, Δt 10 seconds, max(vel): (1.82e-01, 1.81e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 783.161 ms
[ Info: Time: 15.000 minutes, Iteration 270, Δt 10 seconds, max(vel): (1.82e-01, 1.81e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 793.745 ms
[ Info: Time: 15.833 minutes, Iteration 285, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 794.773 ms
[ Info: Time: 16.667 minutes, Iteration 300, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 809.040 ms
[ Info: Time: 17.500 minutes, Iteration 315, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 794.899 ms
[ Info: Time: 18.333 minutes, Iteration 330, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 791.979 ms
[ Info: Time: 19.167 minutes, Iteration 345, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 815.873 ms
[ Info: Time: 20.000 minutes, Iteration 360, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 803.575 ms
[ Info: Time: 20.833 minutes, Iteration 375, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 793.342 ms
[ Info: Time: 21.667 minutes, Iteration 390, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 801.416 ms
[ Info: Time: 22.500 minutes, Iteration 405, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 789.351 ms
[ Info: Time: 23.333 minutes, Iteration 420, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 783.795 ms
[ Info: Time: 24.167 minutes, Iteration 435, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 813.935 ms
[ Info: Time: 25.000 minutes, Iteration 450, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 816.066 ms
[ Info: Time: 25.833 minutes, Iteration 465, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.99, wtime: 826.817 ms
[ Info: Time: 26.667 minutes, Iteration 480, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.98, wtime: 826.248 ms
[ Info: Time: 27.500 minutes, Iteration 495, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.98, wtime: 833.111 ms
[ Info: Time: 28.333 minutes, Iteration 510, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.98, wtime: 833.443 ms
[ Info: Time: 29.167 minutes, Iteration 525, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.98, wtime: 808.598 ms
[ Info: Time: 30.000 minutes, Iteration 540, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.98, wtime: 794.113 ms
[ Info: Time: 30.833 minutes, Iteration 555, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.98, wtime: 811.225 ms
[ Info: Time: 31.667 minutes, Iteration 570, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.98, wtime: 804.547 ms
[ Info: Time: 32.500 minutes, Iteration 585, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.98, wtime: 808.382 ms
[ Info: Time: 33.333 minutes, Iteration 600, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.98, wtime: 821.364 ms
[ Info: Time: 34.167 minutes, Iteration 615, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.98, wtime: 807.599 ms
[ Info: Time: 35.000 minutes, Iteration 630, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.98, wtime: 824.235 ms
[ Info: Time: 35.833 minutes, Iteration 645, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.98, wtime: 817.080 ms
[ Info: Time: 36.667 minutes, Iteration 660, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.98, wtime: 821.432 ms
[ Info: Time: 37.500 minutes, Iteration 675, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.98, wtime: 811.492 ms
[ Info: Time: 38.333 minutes, Iteration 690, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.98, wtime: 815.107 ms
[ Info: Time: 39.167 minutes, Iteration 705, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.98, wtime: 824.988 ms
[ Info: Time: 40.000 minutes, Iteration 720, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 828.434 ms
[ Info: Time: 40.833 minutes, Iteration 735, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 841.348 ms
[ Info: Time: 41.667 minutes, Iteration 750, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 834.679 ms
[ Info: Time: 42.500 minutes, Iteration 765, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 817.290 ms
[ Info: Time: 43.333 minutes, Iteration 780, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 831.116 ms
[ Info: Time: 44.167 minutes, Iteration 795, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 823.776 ms
[ Info: Time: 45.000 minutes, Iteration 810, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 814.016 ms
[ Info: Time: 45.833 minutes, Iteration 825, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 838.235 ms
[ Info: Time: 46.667 minutes, Iteration 840, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 806.626 ms
[ Info: Time: 47.500 minutes, Iteration 855, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 833.379 ms
[ Info: Time: 48.333 minutes, Iteration 870, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 848.702 ms
[ Info: Time: 49.167 minutes, Iteration 885, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 818.269 ms
[ Info: Time: 50.000 minutes, Iteration 900, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 842.887 ms
[ Info: Time: 50.833 minutes, Iteration 915, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 850.640 ms
[ Info: Time: 51.667 minutes, Iteration 930, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 844.567 ms
[ Info: Time: 52.500 minutes, Iteration 945, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 839.848 ms
[ Info: Time: 53.333 minutes, Iteration 960, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.31, min(ℵ): 0.97, wtime: 841.354 ms
[ Info: Time: 54.167 minutes, Iteration 975, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.96, wtime: 826.119 ms
[ Info: Time: 55.000 minutes, Iteration 990, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.96, wtime: 850.083 ms
[ Info: Time: 55.833 minutes, Iteration 1005, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.96, wtime: 852.515 ms
[ Info: Time: 56.667 minutes, Iteration 1020, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.96, wtime: 849.460 ms
[ Info: Time: 57.500 minutes, Iteration 1035, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.96, wtime: 839.462 ms
[ Info: Time: 58.333 minutes, Iteration 1050, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.96, wtime: 845.807 ms
[ Info: Time: 59.167 minutes, Iteration 1065, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.96, wtime: 867.013 ms
[ Info: Time: 1.000 hours, Iteration 1080, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.96, wtime: 834.717 ms
[ Info: Time: 1.014 hours, Iteration 1095, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.96, wtime: 850.094 ms
[ Info: Time: 1.028 hours, Iteration 1110, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.96, wtime: 853.664 ms
[ Info: Time: 1.042 hours, Iteration 1125, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.96, wtime: 843.034 ms
[ Info: Time: 1.056 hours, Iteration 1140, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.96, wtime: 859.517 ms
[ Info: Time: 1.069 hours, Iteration 1155, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.96, wtime: 864.284 ms
[ Info: Time: 1.083 hours, Iteration 1170, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.96, wtime: 844.465 ms
[ Info: Time: 1.097 hours, Iteration 1185, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.96, wtime: 846.775 ms
[ Info: Time: 1.111 hours, Iteration 1200, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.96, wtime: 855.625 ms
[ Info: Time: 1.125 hours, Iteration 1215, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 857.129 ms
[ Info: Time: 1.139 hours, Iteration 1230, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 870.720 ms
[ Info: Time: 1.153 hours, Iteration 1245, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 880.822 ms
[ Info: Time: 1.167 hours, Iteration 1260, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 870.779 ms
[ Info: Time: 1.181 hours, Iteration 1275, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 859.058 ms
[ Info: Time: 1.194 hours, Iteration 1290, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 862.649 ms
[ Info: Time: 1.208 hours, Iteration 1305, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 894.480 ms
[ Info: Time: 1.222 hours, Iteration 1320, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 870.154 ms
[ Info: Time: 1.236 hours, Iteration 1335, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 872.446 ms
[ Info: Time: 1.250 hours, Iteration 1350, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 871.943 ms
[ Info: Time: 1.264 hours, Iteration 1365, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 889.324 ms
[ Info: Time: 1.278 hours, Iteration 1380, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 889.420 ms
[ Info: Time: 1.292 hours, Iteration 1395, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 933.662 ms
[ Info: Time: 1.306 hours, Iteration 1410, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 888.556 ms
[ Info: Time: 1.319 hours, Iteration 1425, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 888.473 ms
[ Info: Time: 1.333 hours, Iteration 1440, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 891.212 ms
[ Info: Time: 1.347 hours, Iteration 1455, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.95, wtime: 903.512 ms
[ Info: Time: 1.361 hours, Iteration 1470, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.94, wtime: 898.014 ms
[ Info: Time: 1.375 hours, Iteration 1485, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.94, wtime: 881.126 ms
[ Info: Time: 1.389 hours, Iteration 1500, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.94, wtime: 889.321 ms
[ Info: Time: 1.403 hours, Iteration 1515, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.94, wtime: 896.079 ms
[ Info: Time: 1.417 hours, Iteration 1530, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.94, wtime: 890.902 ms
[ Info: Time: 1.431 hours, Iteration 1545, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.94, wtime: 892.726 ms
[ Info: Time: 1.444 hours, Iteration 1560, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.94, wtime: 908.962 ms
[ Info: Time: 1.458 hours, Iteration 1575, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.94, wtime: 890.080 ms
[ Info: Time: 1.472 hours, Iteration 1590, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.94, wtime: 905.222 ms
[ Info: Time: 1.486 hours, Iteration 1605, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.94, wtime: 922.769 ms
[ Info: Time: 1.500 hours, Iteration 1620, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.94, wtime: 902.730 ms
[ Info: Time: 1.514 hours, Iteration 1635, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.32, min(ℵ): 0.94, wtime: 915.929 ms
[ Info: Time: 1.528 hours, Iteration 1650, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.94, wtime: 922.065 ms
[ Info: Time: 1.542 hours, Iteration 1665, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.94, wtime: 907.767 ms
[ Info: Time: 1.556 hours, Iteration 1680, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.94, wtime: 934.765 ms
[ Info: Time: 1.569 hours, Iteration 1695, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.94, wtime: 919.537 ms
[ Info: Time: 1.583 hours, Iteration 1710, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.94, wtime: 952.998 ms
[ Info: Time: 1.597 hours, Iteration 1725, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 938.515 ms
[ Info: Time: 1.611 hours, Iteration 1740, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 951.666 ms
[ Info: Time: 1.625 hours, Iteration 1755, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 951.512 ms
[ Info: Time: 1.639 hours, Iteration 1770, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 916.583 ms
[ Info: Time: 1.653 hours, Iteration 1785, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 922.755 ms
[ Info: Time: 1.667 hours, Iteration 1800, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 923.442 ms
[ Info: Time: 1.681 hours, Iteration 1815, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 923.770 ms
[ Info: Time: 1.694 hours, Iteration 1830, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 928.518 ms
[ Info: Time: 1.708 hours, Iteration 1845, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 944.072 ms
[ Info: Time: 1.722 hours, Iteration 1860, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 947.032 ms
[ Info: Time: 1.736 hours, Iteration 1875, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 939.094 ms
[ Info: Time: 1.750 hours, Iteration 1890, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 921.842 ms
[ Info: Time: 1.764 hours, Iteration 1905, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 928.467 ms
[ Info: Time: 1.778 hours, Iteration 1920, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 933.690 ms
[ Info: Time: 1.792 hours, Iteration 1935, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 936.120 ms
[ Info: Time: 1.806 hours, Iteration 1950, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 936.093 ms
[ Info: Time: 1.819 hours, Iteration 1965, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.93, wtime: 913.217 ms
[ Info: Time: 1.833 hours, Iteration 1980, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.92, wtime: 932.283 ms
[ Info: Time: 1.847 hours, Iteration 1995, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.92, wtime: 949.245 ms
[ Info: Time: 1.861 hours, Iteration 2010, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.92, wtime: 927.143 ms
[ Info: Time: 1.875 hours, Iteration 2025, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.92, wtime: 933.750 ms
[ Info: Time: 1.889 hours, Iteration 2040, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.92, wtime: 926.442 ms
[ Info: Time: 1.903 hours, Iteration 2055, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.92, wtime: 933.653 ms
[ Info: Time: 1.917 hours, Iteration 2070, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.92, wtime: 933.152 ms
[ Info: Time: 1.931 hours, Iteration 2085, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.92, wtime: 948.966 ms
[ Info: Time: 1.944 hours, Iteration 2100, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.92, wtime: 961.911 ms
[ Info: Time: 1.958 hours, Iteration 2115, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.92, wtime: 961.534 ms
[ Info: Time: 1.972 hours, Iteration 2130, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.92, wtime: 972.149 ms
[ Info: Time: 1.986 hours, Iteration 2145, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.92, wtime: 939.450 ms
[ Info: Time: 2.000 hours, Iteration 2160, Δt 10 seconds, max(vel): (1.82e-01, 1.82e-01), max(h): 0.33, min(ℵ): 0.92, wtime: 954.650 ms
[ Info: Simulation is stopping after running for 2.310 minutes.
[ Info: Simulation time 2 hours equals or exceeds stop time 2 hours.

Visualize!

Nt = length(htimeseries)
iter = Observable(1)

hi = @lift(htimeseries[$iter][:, :, 1])
ℵi = @lift(ℵtimeseries[$iter][:, :, 1])
ui = @lift(utimeseries[$iter][:, :, 1])
vi = @lift(vtimeseries[$iter][:, :, 1])

fig = Figure()
ax = Axis(fig[1, 1], title = "sea ice thickness")
heatmap!(ax, hi, colormap = :magma,         colorrange = (0.23, 0.37))

ax = Axis(fig[1, 2], title = "sea ice concentration")
heatmap!(ax, ℵi, colormap = Reverse(:deep), colorrange = (0.75, 1))

ax = Axis(fig[2, 1], title = "zonal velocity")
heatmap!(ax, ui, colorrange = (-0.1, 0.1))

ax = Axis(fig[2, 2], title = "meridional velocity")
heatmap!(ax, vi, colorrange = (-0.1, 0.1))

CairoMakie.record(fig, "sea_ice_dynamics.mp4", 1:Nt, framerate = 8) do i
    iter[] = i
    @info "doing iter $i"
end

"sea_ice_dynamics.mp4"

This page was generated using Literate.jl.