From 13399aff622043b6e1b170f4263ad239cbc694f0 Mon Sep 17 00:00:00 2001 From: Simone Silvestri Date: Thu, 23 Jan 2025 16:33:35 +0100 Subject: [PATCH] correction --- ...topologically_conditional_interpolation.jl | 66 ++++++++++++------- 1 file changed, 43 insertions(+), 23 deletions(-) diff --git a/src/Advection/topologically_conditional_interpolation.jl b/src/Advection/topologically_conditional_interpolation.jl index 9c623c6b57..d053afd0a4 100644 --- a/src/Advection/topologically_conditional_interpolation.jl +++ b/src/Advection/topologically_conditional_interpolation.jl @@ -10,18 +10,18 @@ ##### close to the boundary, or a second-order interpolation if i is close to a boundary. ##### -using Oceananigans.Grids: AbstractUnderlyingGrid, Bounded +using Oceananigans.Grids: AbstractUnderlyingGrid, Bounded, RightConnected, LeftConnected, BoundedTopology const AUG = AbstractUnderlyingGrid # Bounded underlying Grids -const AUGX = AUG{<:Any, <:Bounded} -const AUGY = AUG{<:Any, <:Any, <:Bounded} -const AUGZ = AUG{<:Any, <:Any, <:Any, <:Bounded} -const AUGXY = AUG{<:Any, <:Bounded, <:Bounded} -const AUGXZ = AUG{<:Any, <:Bounded, <:Any, <:Bounded} -const AUGYZ = AUG{<:Any, <:Any, <:Bounded, <:Bounded} -const AUGXYZ = AUG{<:Any, <:Bounded, <:Bounded, <:Bounded} +const AUGX = AUG{<:Any, <:BoundedTopology} +const AUGY = AUG{<:Any, <:Any, <:BoundedTopology} +const AUGZ = AUG{<:Any, <:Any, <:Any, <:BoundedTopology} +const AUGXY = AUG{<:Any, <:BoundedTopology, <:BoundedTopology} +const AUGXZ = AUG{<:Any, <:BoundedTopology, <:Any, <:BoundedTopology} +const AUGYZ = AUG{<:Any, <:Any, <:BoundedTopology, <:BoundedTopology} +const AUGXYZ = AUG{<:Any, <:BoundedTopology, <:BoundedTopology, <:BoundedTopology} # Left-biased buffers are smaller by one grid point on the right side; vice versa for right-biased buffers # Center interpolation stencil look at i + 1 (i.e., require one less point on the left) @@ -29,20 +29,40 @@ const AUGXYZ = AUG{<:Any, <:Bounded, <:Bounded, <:Bounded} for dir in (:x, :y, :z) outside_symmetric_haloᶠ = Symbol(:outside_symmetric_halo_, dir, :ᶠ) outside_symmetric_haloᶜ = Symbol(:outside_symmetric_halo_, dir, :ᶜ) - outside_biased_haloᶠ = Symbol(:outside_biased_halo_, dir, :ᶠ) - outside_biased_haloᶜ = Symbol(:outside_biased_halo_, dir, :ᶜ) - required_halo_size = Symbol(:required_halo_size_, dir) + outside_biased_haloᶠ = Symbol(:outside_biased_halo_, dir, :ᶠ) + outside_biased_haloᶜ = Symbol(:outside_biased_halo_, dir, :ᶜ) + required_halo_size = Symbol(:required_halo_size_, dir) @eval begin - @inline $outside_symmetric_haloᶠ(i, N, adv) = (i >= $required_halo_size(adv) + 1) & (i <= N + 1 - $required_halo_size(adv)) - @inline $outside_symmetric_haloᶜ(i, N, adv) = (i >= $required_halo_size(adv)) & (i <= N + 1 - $required_halo_size(adv)) - - @inline $outside_biased_haloᶠ(i, N, adv) = (i >= $required_halo_size(adv) + 1) & (i <= N + 1 - ($required_halo_size(adv) - 1)) & # Left bias - (i >= $required_halo_size(adv)) & (i <= N + 1 - $required_halo_size(adv)) # Right bias - @inline $outside_biased_haloᶜ(i, N, adv) = (i >= $required_halo_size(adv)) & (i <= N + 1 - ($required_halo_size(adv) - 1)) & # Left bias - (i >= $required_halo_size(adv) - 1) & (i <= N + 1 - $required_halo_size(adv)) # Right bias + # Bounded topologies + @inline $outside_symmetric_haloᶠ(i, ::Type{Bounded}, N, adv) = (i >= $required_halo_size(adv) + 1) & (i <= N + 1 - $required_halo_size(adv)) + @inline $outside_symmetric_haloᶜ(i, ::Type{Bounded}, N, adv) = (i >= $required_halo_size(adv)) & (i <= N + 1 - $required_halo_size(adv)) + + @inline $outside_biased_haloᶠ(i, ::Type{Bounded}, N, adv) = (i >= $required_halo_size(adv) + 1) & (i <= N + 1 - ($required_halo_size(adv) - 1)) & # Left bias + (i >= $required_halo_size(adv)) & (i <= N + 1 - $required_halo_size(adv)) # Right bias + @inline $outside_biased_haloᶜ(i, ::Type{Bounded}, N, adv) = (i >= $required_halo_size(adv)) & (i <= N + 1 - ($required_halo_size(adv) - 1)) & # Left bias + (i >= $required_halo_size(adv) - 1) & (i <= N + 1 - $required_halo_size(adv)) # Right bias + + # Right connected topologies (only test the left side, i.e. the bounded side) + @inline $outside_symmetric_haloᶠ(i, ::Type{RightConnected}, N, adv) = i >= $required_halo_size(adv) + 1 + @inline $outside_symmetric_haloᶜ(i, ::Type{RightConnected}, N, adv) = i >= $required_halo_size(adv) + + @inline $outside_biased_haloᶠ(i, ::Type{RightConnected}, N, adv) = (i >= $required_halo_size(adv) + 1) & # Left bias + (i >= $required_halo_size(adv)) # Right bias + @inline $outside_biased_haloᶜ(i, ::Type{RightConnected}, N, adv) = (i >= $required_halo_size(adv)) & # Left bias + (i >= $required_halo_size(adv) - 1) # Right bias + + # Left bounded topologies (only test the right side, i.e. the bounded side) + @inline $outside_symmetric_haloᶠ(i, ::Type{LeftConnected}, N, adv) = (i <= N + 1 - $required_halo_size(adv)) + @inline $outside_symmetric_haloᶜ(i, ::Type{LeftConnected}, N, adv) = (i <= N + 1 - $required_halo_size(adv)) + + @inline $outside_biased_haloᶠ(i, ::Type{LeftConnected}, N, adv) = (i <= N + 1 - ($required_halo_size(adv) - 1)) & # Left bias + (i <= N + 1 - $required_halo_size(adv)) # Right bias + @inline $outside_biased_haloᶜ(i, ::Type{LeftConnected}, N, adv) = (i <= N + 1 - ($required_halo_size(adv) - 1)) & # Left bias + (i <= N + 1 - $required_halo_size(adv)) # Right bias end end + # Separate High order advection from low order advection const HOADV = Union{WENO, Tuple(Centered{N} for N in advection_buffers[2:end])..., @@ -76,21 +96,21 @@ for bias in (:symmetric, :biased) if ξ == :x @eval begin @inline $alt1_interp(i, j, k, grid::AUGX, scheme::HOADV, args...) = - ifelse($outside_buffer(i, grid.Nx, scheme), + ifelse($outside_buffer(i, topology(grid, 1), grid.Nx, scheme), $interp(i, j, k, grid, scheme, args...), $alt2_interp(i, j, k, grid, scheme.buffer_scheme, args...)) end elseif ξ == :y @eval begin @inline $alt1_interp(i, j, k, grid::AUGY, scheme::HOADV, args...) = - ifelse($outside_buffer(j, grid.Ny, scheme), + ifelse($outside_buffer(j, topology(grid, 2), grid.Ny, scheme), $interp(i, j, k, grid, scheme, args...), $alt2_interp(i, j, k, grid, scheme.buffer_scheme, args...)) end elseif ξ == :z @eval begin @inline $alt1_interp(i, j, k, grid::AUGZ, scheme::HOADV, args...) = - ifelse($outside_buffer(k, grid.Nz, scheme), + ifelse($outside_buffer(k, topology(grid, 3), grid.Nz, scheme), $interp(i, j, k, grid, scheme, args...), $alt2_interp(i, j, k, grid, scheme.buffer_scheme, args...)) end @@ -100,11 +120,11 @@ for bias in (:symmetric, :biased) end @inline _multi_dimensional_reconstruction_x(i, j, k, grid::AUGX, scheme, interp, args...) = - ifelse(outside_symmetric_bufferᶜ(i, grid.Nx, scheme), + ifelse(outside_symmetric_bufferᶜ(i, topology(grid, 1), grid.Nx, scheme), multi_dimensional_reconstruction_x(i, j, k, grid::AUGX, scheme, interp, args...), interp(i, j, k, grid, scheme, args...)) @inline _multi_dimensional_reconstruction_y(i, j, k, grid::AUGY, scheme, interp, args...) = - ifelse(outside_symmetric_bufferᶜ(j, grid.Ny, scheme), + ifelse(outside_symmetric_bufferᶜ(j, topology(grid, 2), grid.Ny, scheme), multi_dimensional_reconstruction_y(i, j, k, grid::AUGY, scheme, interp, args...), interp(i, j, k, grid, scheme, args...))