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.

Sign up for GitHub

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

Jump to bottom

Split-Explicit free-surface solver interface revisit #3311

Closed
wants to merge 11 commits into from
Closed

Split-Explicit free-surface solver interface revisit #3311

Changes from 1 commit
Commits
Show all changes
11 commits
Select commit Hold shift + click to select a range
afb4661
bit of cleanup
navidcy Oct 5, 2023
ee393b4
Merge branch 'main' into ncc-glw-ss/split-explicit-interface
navidcy Dec 29, 2023
c876dff
Merge branch 'main' into ncc-glw-ss/split-explicit-interface
navidcy Jan 16, 2024
44e10c0
use grid as arg in FixedTimeStepSize constructor
navidcy Jan 16, 2024
e545fe5
code alignment
navidcy Jan 16, 2024
3e30dd8
remove stray s
navidcy Jan 16, 2024
cd11a69
fix show method for SplitExplicitFreeSurface
navidcy Jan 16, 2024
edca7a6
fix typos
navidcy Jan 16, 2024
2faa01f
beautifications
navidcy Jan 16, 2024
aa2c8cd
simplify
navidcy Jan 16, 2024
b7ef2c4
remove stray spaces
navidcy Jan 16, 2024
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
25 changes: 14 additions & 11 deletions src/Models/HydrostaticFreeSurfaceModels/split_explicit_free_surface.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -47,8 +47,8 @@ Keyword Arguments

- `cfl`: If set then the number of `substeps` are computed based on the advective timescale imposed from the
barotropic gravity-wave speed, computed with depth `grid.Lz`. If `fixed_Δt` is provided then the number of
`substeps` will adapt to maintain an exact cfl. If not the effective cfl will be always lower than the
specified `cfl` provided that the baroclinic time step `Δt_baroclinic < fixed_Δt`
`substeps` will adapt to maintain an exact cfl. If not, the effective cfl will always be lower than the
specified `cfl` provided that the baroclinic time step `Δt_baroclinic < fixed_Δt`.

!!! info "Needed keyword arguments"
Either `substeps` _or_ `cfl` (with `grid`) need to be prescribed.
Expand Down Expand Up @@ -91,12 +91,14 @@ function SplitExplicitFreeSurface(grid; gravitational_acceleration = g_Earth,
if eltype(settings) != eltype(grid)
@warn "Using $(eltype(settings)) settings for the SplitExplicitFreeSurface on a $(eltype(grid)) grid"
end

η = ZFaceField(grid, indices = (:, :, size(grid, 3)+1))

gravitational_acceleration = convert(eltype(grid), gravitational_acceleration)

return SplitExplicitFreeSurface(η, SplitExplicitState(grid), SplitExplicitAuxiliaryFields(grid),
gravitational_acceleration, settings)
return SplitExplicitFreeSurface(η, SplitExplicitState(grid),
SplitExplicitAuxiliaryFields(grid),
gravitational_acceleration, settings)
end

"""
Expand Down Expand Up @@ -206,7 +208,7 @@ function SplitExplicitAuxiliaryFields(grid::AbstractGrid)

dz = GridMetricOperation((Face, Center, Center), Δz, grid)
sum!(Hᶠᶜ, dz)

dz = GridMetricOperation((Center, Face, Center), Δz, grid)
sum!(Hᶜᶠ, dz)

Expand All @@ -216,7 +218,7 @@ function SplitExplicitAuxiliaryFields(grid::AbstractGrid)
fill_halo_regions!((Hᶠᶜ, Hᶜᶠ, Hᶜᶜ))

kernel_parameters = :xy

return SplitExplicitAuxiliaryFields(Gᵁ, Gⱽ, Hᶠᶜ, Hᶜᶠ, Hᶜᶜ, kernel_parameters)
end

Expand Down Expand Up @@ -261,15 +263,16 @@ end

function FixedTimeStepSize(FT::DataType = Float64;
cfl = 0.7,
grid,
grid,
averaging_kernel = averaging_shape_function,
gravitational_acceleration = g_Earth)

Δx⁻² = topology(grid)[1] == Flat ? 0 : 1 / minimum_xspacing(grid)^2
Δy⁻² = topology(grid)[2] == Flat ? 0 : 1 / minimum_yspacing(grid)^2
Δs = sqrt(1 / (Δx⁻² + Δy⁻²))
Comment on lines 269 to 271
Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@simone-silvestri is this the proper way to compute the spacing for the CFL?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think we should make this a small function on its own since users might wanna use this in their scripts?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

the function should handle Flat ?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

minimum_xspacing should handle Flat. You shouldn't need to check manually.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

True!!

Copy link
Collaborator Author

@navidcy navidcy Oct 5, 2023 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Hm. I think we stumbled on this before. What’s the miminum xspacing for Flat x? 0 I presume.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Inf

Copy link
Member

@glwagner glwagner Oct 5, 2023 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

this use case illustrates it perfectly. The neutral value for a minimum calculation is Inf (everything is less than Inf).

Copy link
Collaborator Author

@navidcy navidcy Jan 16, 2024 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

At the moment, for Flat dimensions the spacings have the value f 1 in the grids:

Oceananigans.jl/src/Grids/grid_generation.jl

Lines 124 to 126 in 1c2a6f8

# Flat domains
generate_coordinate(FT, ::Flat, N, H, coord::Tuple{<:Number, <:Number}, coordinate_name, arch) =
FT(1), range(1, 1, length=N), range(1, 1, length=N), FT(1), FT(1)

Then, when x/y/zminimum_spacing is called for those Flat dimension the 1 is returned...

We should really have 0 or Inf as the spacing for Flat dimensions, right?


wave_speed = sqrt(gravitational_acceleration * grid.Lz)

Δt_barotropic = convert(FT, cfl * Δs / wave_speed)

return FixedTimeStepSize(Δt_barotropic, averaging_kernel)
Expand All @@ -284,7 +287,7 @@ end
idx = searchsortedlast(averaging_weights, 0, rev=true)
substeps = idx

averaging_weights = averaging_weights[1:idx]
averaging_weights = @inbounds averaging_weights[1:idx]
averaging_weights ./= sum(averaging_weights)

return Δτ, averaging_weights
Expand All @@ -298,7 +301,7 @@ function SplitExplicitSettings(FT::DataType=Float64;
gravitational_acceleration = g_Earth,
averaging_kernel = averaging_shape_function,
timestepper = ForwardBackwardScheme())

if (!isnothing(substeps) && !isnothing(cfl)) || (isnothing(substeps) && isnothing(cfl))
throw(ArgumentError("either specify a cfl or a number of substeps"))
end
Expand Down