-
Notifications
You must be signed in to change notification settings - Fork 194
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.
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
Fix _immered_cell
and adapt for PartialCellBottom
#3682
Conversation
How does this PR change |
# Face node below current cell | ||
z = znode(i, j, k, underlying_grid, c, c, f) | ||
h = @inbounds ib.bottom_height[i, j, 1] | ||
ϵ = ib.minimum_fractional_cell_height | ||
# z + Δz is equal to the face above the current cell | ||
Δz = Δzᶜᶜᶜ(i, j, k, ibg.underlying_grid) | ||
return z + Δz * ϵ ≤ h |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Should we expect this change to fix the weird behaviors we were seeing in #3529 ?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I am not sure, @jm-c is the author of this fix, and he said that with this fix the immersed boundary example did not NaN immediately but it still NaNed after a little bit so we believe there is more bugs to be squashed. You can try that example on this and see if it makes any difference. It will most likely change the solution.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I am happy to try the example from #3529 it on this branch, when you think it's ready for testing.
What are we waiting for to merge this PR? |
I guess to have the partial cells working properly, @jm-c was saying there is still some fixes to be done. |
The first thing maybe is to add a test that catches the bug. |
@ali-ramadhan does this PR fix your issue? I will try the internal tide example on this PR, add some tests for |
This is the status of the internal_tide example on this PR, not producing NaNs anymore, but still not really working. internal_tide.mp4 |
@simone-silvestri I think another error has now come up with the original MWE in issue #3681. I tried the MWE: using Oceananigans
using Oceananigans.ImmersedBoundaries: PartialCellBottom
underlying_grid = RectilinearGrid(GPU(), topology=(Periodic, Flat, Bounded), size=(4, 4), x=(0, 1), z=(0, 1))
slope(x) = x
grid = ImmersedBoundaryGrid(underlying_grid, PartialCellBottom(slope))
model = HydrostaticFreeSurfaceModel(; grid) Error: ERROR: CUDA error: an illegal memory access was encountered (code 700, ERROR_ILLEGAL_ADDRESS)
Stacktrace:
[1] throw_api_error(res::CUDA.cudaError_enum)
@ CUDA ~/.julia/packages/CUDA/Tl08O/lib/cudadrv/libcuda.jl:30
[2] check
@ ~/.julia/packages/CUDA/Tl08O/lib/cudadrv/libcuda.jl:37 [inlined]
[3] cuStreamGetCaptureInfo
@ ~/.julia/packages/CUDA/Tl08O/lib/utils/call.jl:34 [inlined]
[4] capture_status(stream::CUDA.CuStream)
@ CUDA ~/.julia/packages/CUDA/Tl08O/lib/cudadrv/graph.jl:174
[5] is_capturing
@ ~/.julia/packages/CUDA/Tl08O/lib/cudadrv/graph.jl:179 [inlined]
[6] convert(::Type{CUDA.CuPtr{Float64}}, managed::CUDA.Managed{CUDA.DeviceMemory})
@ CUDA ~/.julia/packages/CUDA/Tl08O/src/memory.jl:539
[7] unsafe_convert
@ ~/.julia/packages/CUDA/Tl08O/src/array.jl:434 [inlined]
[8] #pointer#1123
@ ~/.julia/packages/CUDA/Tl08O/src/array.jl:392 [inlined]
[9] pointer (repeats 2 times)
@ ~/.julia/packages/CUDA/Tl08O/src/array.jl:384 [inlined]
[10] unsafe_convert(::Type{CUDA.CuDeviceArray{Float64, 3, 1}}, a::CUDA.CuArray{Float64, 3, CUDA.DeviceMemory})
@ CUDA ~/.julia/packages/CUDA/Tl08O/src/array.jl:454
[11] adapt_storage
@ ~/.julia/packages/CUDA/Tl08O/src/compiler/execution.jl:162 [inlined]
[12] adapt_structure
@ ~/.julia/packages/Adapt/7T9au/src/Adapt.jl:57 [inlined]
[13] adapt
@ ~/.julia/packages/Adapt/7T9au/src/Adapt.jl:40 [inlined]
[14] #1
@ ~/.julia/packages/OffsetArrays/hwmnB/ext/OffsetArraysAdaptExt.jl:9 [inlined]
[15] parent_call
@ ~/.julia/packages/OffsetArrays/hwmnB/src/OffsetArrays.jl:315 [inlined]
[16] adapt_structure
@ ~/.julia/packages/OffsetArrays/hwmnB/ext/OffsetArraysAdaptExt.jl:9 [inlined]
[17] adapt
@ ~/.julia/packages/Adapt/7T9au/src/Adapt.jl:40 [inlined]
[18] cudaconvert
@ ~/.julia/packages/CUDA/Tl08O/src/compiler/execution.jl:198 [inlined]
[19] map
@ ./tuple.jl:294 [inlined]
[20] map(f::typeof(CUDA.cudaconvert), t::Tuple{KernelAbstractions.CompilerMetadata{…}, OffsetArrays.OffsetArray{…}, Nothing, Nothing, Tuple{…}, ImmersedBoundaryGrid{…}, Tuple{}})
@ Base ./tuple.jl:294
[21] macro expansion
@ ~/.julia/packages/CUDA/Tl08O/src/compiler/execution.jl:110 [inlined]
[22] (::KernelAbstractions.Kernel{…})(::OffsetArrays.OffsetArray{…}, ::Vararg{…}; ndrange::Nothing, workgroupsize::Nothing)
@ CUDA.CUDAKernels ~/.julia/packages/CUDA/Tl08O/src/CUDAKernels.jl:103
[23] (::KernelAbstractions.Kernel{CUDA.CUDAKernels.CUDABackend, KernelAbstractions.NDIteration.StaticSize{…}, KernelAbstractions.NDIteration.StaticSize{…}, typeof(Oceananigans.BoundaryConditions.gpu__no_fill!)})(::OffsetArrays.OffsetArray{Float64, 3, CUDA.CuArray{…}}, ::Vararg{Any})
@ CUDA.CUDAKernels ~/.julia/packages/CUDA/Tl08O/src/CUDAKernels.jl:89
[24] launch!(::GPU, ::ImmersedBoundaryGrid{…}, ::Tuple{…}, ::typeof(Oceananigans.BoundaryConditions._no_fill!), ::OffsetArrays.OffsetArray{…}, ::Vararg{…}; include_right_boundaries::Bool, reduced_dimensions::Tuple{}, location::Nothing, active_cells_map::Nothing, kwargs::@Kwargs{})
@ Oceananigans.Utils /global/u1/a/alir/Oceananigans.jl/src/Utils/kernel_launching.jl:168
[25] launch!
@ /global/u1/a/alir/Oceananigans.jl/src/Utils/kernel_launching.jl:154 [inlined]
[26] #fill_open_boundary_regions!#48
@ /global/u1/a/alir/Oceananigans.jl/src/BoundaryConditions/fill_halo_regions_open.jl:20 [inlined]
[27] fill_open_boundary_regions!
@ /global/u1/a/alir/Oceananigans.jl/src/BoundaryConditions/fill_halo_regions_open.jl:9 [inlined]
[28] fill_halo_regions!(::OffsetArrays.OffsetArray{…}, ::FieldBoundaryConditions{…}, ::Tuple{…}, ::Tuple{…}, ::ImmersedBoundaryGrid{…}; fill_boundary_normal_velocities::Bool, kwargs::@Kwargs{…})
@ Oceananigans.BoundaryConditions /global/u1/a/alir/Oceananigans.jl/src/BoundaryConditions/fill_halo_regions.jl:53
[29] fill_halo_regions!
@ /global/u1/a/alir/Oceananigans.jl/src/BoundaryConditions/fill_halo_regions.jl:48 [inlined]
[30] #fill_halo_regions!#61
@ /global/u1/a/alir/Oceananigans.jl/src/Fields/field.jl:773 [inlined]
[31] fill_halo_regions!(::Field{Face, Center, Nothing, Nothing, ImmersedBoundaryGrid{…}, Tuple{…}, OffsetArrays.OffsetArray{…}, Float64, FieldBoundaryConditions{…}, Nothing, Oceananigans.Fields.FieldBoundaryBuffers{…}})
@ Oceananigans.Fields /global/u1/a/alir/Oceananigans.jl/src/Fields/field.jl:770
[32] fill_halo_regions!(::Tuple{Field{…}, Field{…}}; kwargs::@Kwargs{})
@ Oceananigans.Fields /global/u1/a/alir/Oceananigans.jl/src/Fields/field_tuples.jl:64
[33] fill_halo_regions!(::Tuple{Field{…}, Field{…}})
@ Oceananigans.Fields /global/u1/a/alir/Oceananigans.jl/src/Fields/field_tuples.jl:56
[34] Oceananigans.Models.HydrostaticFreeSurfaceModels.SplitExplicitAuxiliaryFields(grid::ImmersedBoundaryGrid{Float64, Periodic, Flat, Bounded, RectilinearGrid{…}, PartialCellBottom{…}, Nothing, Nothing, GPU})
@ Oceananigans.Models.HydrostaticFreeSurfaceModels /global/u1/a/alir/Oceananigans.jl/src/Models/HydrostaticFreeSurfaceModels/split_explicit_free_surface.jl:239
[35] materialize_free_surface(free_surface::SplitExplicitFreeSurface{…}, velocities::@NamedTuple{…}, grid::ImmersedBoundaryGrid{…})
@ Oceananigans.Models.HydrostaticFreeSurfaceModels /global/u1/a/alir/Oceananigans.jl/src/Models/HydrostaticFreeSurfaceModels/split_explicit_free_surface.jl:121
[36] HydrostaticFreeSurfaceModel(; grid::ImmersedBoundaryGrid{…}, clock::Clock{…}, momentum_advection::Centered{…}, tracer_advection::Centered{…}, buoyancy::Nothing, coriolis::Nothing, free_surface::SplitExplicitFreeSurface{…}, tracers::Nothing, forcing::@NamedTuple{}, closure::Nothing, boundary_conditions::@NamedTuple{}, particles::Nothing, biogeochemistry::Nothing, velocities::Nothing, pressure::Nothing, diffusivity_fields::Nothing, auxiliary_fields::@NamedTuple{})
@ Oceananigans.Models.HydrostaticFreeSurfaceModels /global/u1/a/alir/Oceananigans.jl/src/Models/HydrostaticFreeSurfaceModels/hydrostatic_free_surface_model.jl:177
[37] top-level scope
@ REPL[6]:1
Some type information was truncated. Use `show(err)` to see complete types. |
Hmmm, I cannot reproduce it on the CPU with Oceananigans.jl/src/ImmersedBoundaries/partial_cell_bottom.jl Lines 155 to 161 in 0fd3986
|
I was able to reproduce on the CPU with Stacktrace from
|
Confirming that the MWE just above works on this branch on both CPU and GPU! |
Also confirming that my original MWE from #3681 also works on this branch now! using Oceananigans
using Oceananigans.ImmersedBoundaries: PartialCellBottom
underlying_grid = RectilinearGrid(GPU(), topology=(Periodic, Flat, Bounded), size=(4, 4), x=(0, 1), z=(0, 1))
slope(x) = x
grid = ImmersedBoundaryGrid(underlying_grid, PartialCellBottom(slope))
model = HydrostaticFreeSurfaceModel(; grid)
time_step!(model, 1) |
What's keeping us from merging this? It fixes something so its a step forward, I think we should merge it. |
…/adapt-bottom-height
I implemented the clamping from #3787 and it seems to have fixed the issue: GridFittedBottominternal_tide_gf.mp4PartialCellBottominternal_tide_pc.mp4I'm not sure I really understand why it would make a difference. I'm wondering @simone-silvestri are you sure that your test was using this code? Either way I see no reason not to merge this. Perhaps we will continue to find improvements for PartialCellBottom but I think this is ready to be used in an example. Note that in the internal tide case the partial cells have the interesting effect of weakening the vertical velocity. That could make sense if the topography is somehow "less effectively steep" than in the full cell case. |
This PR now supercedes #3787 |
Now we have merged, should I revisit #3529? |
It's a good idea! Also @apaloczy had some tests that might be good to try to re-run. |
PartialCellBottom
_immered_cell
and adapt for PartialCellBottom
I suspect it won't change anything though because this PR fixes a bug that affects when criteria about the minimum fractional cell height is met. I don't think that criteria is met in your example (which is good) |
closes #3681
We can probably use this PR to add also a validation for the
PartialCellBottom
and the bug-fixes in the implementation.What do people think?
cc @jm-c