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.

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

Add test for autodifferentiating hydrostatic turbulence simulation #3867

Open
wants to merge 37 commits into
base: main
Choose a base branch
from

Conversation

glwagner
Copy link
Member

@glwagner glwagner commented Oct 24, 2024

This PR adds a test that automatically differentiates a hydrostatic turbulence simulation that uses Centered momentum advection and an ExplicitFreeSurface.

@wsmoses @jlk9

Previously we had deduced that this model was differentiable, but it looks like there are still some issues. Not sure what the other test was doing. I think we should consider this separately from #3822 since that PR is currently working on the SplitExplicitFreeSurface.

@glwagner
Copy link
Member Author

Here's the error that I get right now when running this locally:

ERROR: LoadError: Enzyme cannot deduce type
Current scope:
; Function Attrs: mustprogress nofree readonly willreturn
define "enzyme_type"="{[0]:Integer, [1]:Integer, [2]:Integer, [3]:Integer, [4]:Integer, [5]:Integer, [6]:Integer, [7]:Integer, [8]:Integer, [9]:Integer, [10]:Integer, [11]:Integer, [12]:Integer, [13]:Integer, [14]:Integer, [15]:Integer, [16]:Integer, [17]:Integer, [18]:Integer, [19]:Integer, [20]:Integer, [21]:Integer, [22]:Integer, [23]:Integer, [24]:Integer, [25]:Integer, [26]:Integer, [27]:Integer, [28]:Integer, [29]:Integer, [30]:Integer, [31]:Integer, [32]:Integer, [33]:Integer, [34]:Integer, [35]:Integer, [36]:Integer, [37]:Integer, [38]:Integer, [39]:Integer, [40]:Integer, [41]:Integer, [42]:Integer, [43]:Integer, [44]:Integer, [45]:Integer, [46]:Integer, [47]:Integer, [48]:Float@double, [56]:Float@double, [64]:Float@double, [72]:Float@double, [80]:Float@double, [88]:Float@double, [96]:Float@double, [104]:Float@double, [112]:Float@double, [120]:Integer, [121]:Integer, [122]:Integer, [123]:Integer, [124]:Integer, [125]:Integer, [126]:Integer, [127]:Integer, [128]:Integer, [129]:Integer, [130]:Integer, [131]:Integer, [132]:Integer, [133]:Integer, [134]:Integer, [135]:Integer, [136]:Integer, [137]:Integer, [138]:Integer, [139]:Integer, [140]:Integer, [141]:Integer, [142]:Integer, [143]:Integer, [144]:Float@double, [152]:Float@double, [160]:Float@double, [168]:Float@double, [176]:Integer, [177]:Integer, [178]:Integer, [179]:Integer, [180]:Integer, [181]:Integer, [182]:Integer, [183]:Integer, [184]:Integer, [185]:Integer, [186]:Integer, [187]:Integer, [188]:Integer, [189]:Integer, [190]:Integer, [191]:Integer, [192]:Integer, [193]:Integer, [194]:Integer, [195]:Integer, [196]:Integer, [197]:Integer, [198]:Integer, [199]:Integer, [200]:Float@double, [208]:Float@double, [216]:Float@double, [224]:Float@double, [232]:Float@double, [240]:Float@double, [248]:Integer, [249]:Integer, [250]:Integer, [251]:Integer, [252]:Integer, [253]:Integer, [254]:Integer, [255]:Integer, [256]:Integer, [257]:Integer, [258]:Integer, [259]:Integer, [260]:Integer, [261]:Integer, [262]:Integer, [263]:Integer, [264]:Integer, [265]:Integer, [266]:Integer, [267]:Integer, [268]:Integer, [269]:Integer, [270]:Integer, [271]:Integer, [272]:Float@double, [280]:Float@double, [288]:Float@double, [296]:Float@double, [304]:Integer, [305]:Integer, [306]:Integer, [307]:Integer, [308]:Integer, [309]:Integer, [310]:Integer, [311]:Integer, [312]:Integer, [313]:Integer, [314]:Integer, [315]:Integer, [316]:Integer, [317]:Integer, [318]:Integer, [319]:Integer, [320]:Integer, [321]:Integer, [322]:Integer, [323]:Integer, [324]:Integer, [325]:Integer, [326]:Integer, [327]:Integer, [328]:Float@double, [336]:Float@double, [344]:Float@double, [352]:Float@double, [360]:Float@double, [368]:Float@double, [376]:Integer, [377]:Integer, [378]:Integer, [379]:Integer, [380]:Integer, [381]:Integer, [382]:Integer, [383]:Integer, [384]:Integer, [385]:Integer, [386]:Integer, [387]:Integer, [388]:Integer, [389]:Integer, [390]:Integer, [391]:Integer, [392]:Integer, [393]:Integer, [394]:Integer, [395]:Integer, [396]:Integer, [397]:Integer, [398]:Integer, [399]:Integer, [400]:Float@double, [408]:Float@double, [416]:Float@double, [424]:Float@double, [432]:Integer, [433]:Integer, [434]:Integer, [435]:Integer, [436]:Integer, [437]:Integer, [438]:Integer, [439]:Integer, [440]:Integer, [441]:Integer, [442]:Integer, [443]:Integer, [444]:Integer, [445]:Integer, [446]:Integer, [447]:Integer, [448]:Integer, [449]:Integer, [450]:Integer, [451]:Integer, [452]:Integer, [453]:Integer, [454]:Integer, [455]:Integer, [456]:Pointer, [456,0]:Pointer, [456,0,-1]:Float@double, [456,8]:Integer, [456,9]:Integer, [456,10]:Integer, [456,11]:Integer, [456,12]:Integer, [456,13]:Integer, [456,14]:Integer, [456,15]:Integer, [456,16]:Integer, [456,17]:Integer, [456,18]:Integer, [456,19]:Integer, [456,20]:Integer, [456,21]:Integer, [456,22]:Integer, [456,23]:Integer, [456,24]:Integer, [456,25]:Integer, [456,26]:Integer, [456,27]:Integer, [456,28]:Integer, [456,29]:Integer, [456,30]:Integer, [456,31]:Integer, [456,32]:Integer, [456,33]:Integer, [456,34]:Integer, [456,35]:Integer, [456,36]:Integer, [456,37]:Integer, [456,38]:Integer, [456,39]:Integer, [464]:Integer, [465]:Integer, [466]:Integer, [467]:Integer, [468]:Integer, [469]:Integer, [470]:Integer, [471]:Integer, [472]:Integer, [473]:Integer, [474]:Integer, [475]:Integer, [476]:Integer, [477]:Integer, [478]:Integer, [479]:Integer, [480]:Integer, [481]:Integer, [482]:Integer, [483]:Integer, [484]:Integer, [485]:Integer, [486]:Integer, [487]:Integer, [488]:Pointer, [496]:Integer, [497]:Integer, [498]:Integer, [499]:Integer, [500]:Integer}" "enzymejl_parmtype"="5997836240" "enzymejl_parmtype_ref"="1" [2 x { { i64, i64, i64, i64, i64, i64, double, double, double, double, double, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, double, double, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, double, double, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] } }, { {} addrspace(10)*, [3 x i64] }, {} addrspace(10)* }] @preprocess_julia_filter_21024_inner.1({ { { i64, i64, i64, i64, i64, i64, double, double, double, double, double, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, double, double, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, double, double, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] } }, { {} addrspace(10)*, [3 x i64] }, {} addrspace(10)* }, { { i64, i64, i64, i64, i64, i64, double, double, double, double, double, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, double, double, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, double, double, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] } }, { {} addrspace(10)*, [3 x i64] }, {} addrspace(10)* }, { { i64, i64, i64, i64, i64, i64, double, double, double, double, double, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, double, double, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, double, double, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] } }, { {} addrspace(10)*, [3 x i64] }, {} addrspace(10)*, { [2 x i64] } } } addrspace(11)* noalias nocapture nofree noundef nonnull readonly align 8 dereferenceable(1504) "enzyme_type"="{[-1]:Pointer, [-1,0]:Integer, [-1,1]:Integer, [-1,2]:Integer, [-1,3]:Integer, [-1,4]:Integer, [-1,5]:Integer, [-1,6]:Integer, [-1,7]:Integer, [-1,8]:Integer, [-1,9]:Integer, [-1,10]:Integer, [-1,11]:Integer, [-1,12]:Integer, [-1,13]:Integer, [-1,14]:Integer, [-1,15]:Integer, [-1,16]:Integer, [-1,17]:Integer, [-1,18]:Integer, [-1,19]:Integer, [-1,20]:Integer, [-1,21]:Integer, [-1,22]:Integer, [-1,23]:Integer, [-1,24]:Integer, [-1,25]:Integer, [-1,26]:Integer, [-1,27]:Integer, [-1,28]:Integer, [-1,29]:Integer, [-1,30]:Integer, [-1,31]:Integer, [-1,32]:Integer, [-1,33]:Integer, [-1,34]:Integer, [-1,35]:Integer, [-1,36]:Integer, [-1,37]:Integer, [-1,38]:Integer, [-1,39]:Integer, [-1,40]:Integer, [-1,41]:Integer, [-1,42]:Integer, [-1,43]:Integer, [-1,44]:Integer, [-1,45]:Integer, [-1,46]:Integer, [-1,47]:Integer, [-1,48]:Float@double, [-1,56]:Float@double, [-1,64]:Float@double, [-1,72]:Float@double, [-1,80]:Float@double, [-1,88]:Float@double, [-1,96]:Float@double, [-1,104]:Float@double, [-1,112]:Float@double, [-1,120]:Integer, [-1,121]:Integer, [-1,122]:Integer, [-1,123]:Integer, [-1,124]:Integer, [-1,125]:Integer, [-1,126]:Integer, [-1,127]:Integer, [-1,128]:Integer, [-1,129]:Integer, [-1,130]:Integer, [-1,131]:Integer, [-1,132]:Integer, [-1,133]:Integer, [-1,134]:Integer, [-1,135]:Integer, [-1,136]:Integer, [-1,137]:Integer, [-1,138]:Integer, [-1,139]:Integer, [-1,140]:Integer, [-1,141]:Integer, [-1,142]:Integer, [-1,143]:Integer, [-1,144]:Float@double, [-1,152]:Float@double, [-1,160]:Float@double, [-1,168]:Float@double, [-1,176]:Integer, [-1,177]:Integer, [-1,178]:Integer, [-1,179]:Integer, [-1,180]:Integer, [-1,181]:Integer, [-1,182]:Integer, [-1,183]:Integer, [-1,184]:Integer, [-1,185]:Integer, [-1,186]:Integer, [-1,187]:Integer, [-1,188]:Integer, [-1,189]:Integer, [-1,190]:Integer, [-1,191]:Integer, [-1,192]:Integer, [-1,193]:Integer, [-1,194]:Integer, [-1,195]:Integer, [-1,196]:Integer, [-1,197]:Integer, [-1,198]:Integer, [-1,199]:Integer, [-1,200]:Float@double, [-1,208]:Float@double, [-1,216]:Float@double, [-1,224]:Float@double, [-1,232]:Float@double, [-1,240]:Float@double, [-1,248]:Integer, [-1,249]:Integer, [-1,250]:Integer, [-1,251]:Integer, [-1,252]:Integer, [-1,253]:Integer, [-1,254]:Integer, [-1,255]:Integer, [-1,256]:Integer, [-1,257]:Integer, [-1,258]:Integer, [-1,259]:Integer, [-1,260]:Integer, [-1,261]:Integer, [-1,262]:Integer, [-1,263]:Integer, [-1,264]:Integer, [-1,265]:Integer, [-1,266]:Integer, [-1,267]:Integer, [-1,268]:Integer, [-1,269]:Integer, [-1,270]:Integer, [-1,271]:Integer, [-1,272]:Float@double, [-1,280]:Float@double, [-1,288]:Float@double, [-1,296]:Float@double, [-1,304]:Integer, [-1,305]:Integer, [-1,306]:Integer, [-1,307]:Integer, [-1,308]:Integer, [-1,309]:Integer, [-1,310]:Integer, [-1,311]:Integer, [-1,312]:Integer, [-1,313]:Integer, [-1,314]:Integer, [-1,315]:Integer, [-1,316]:Integer, [-1,317]:Integer, [-1,318]:Integer, [-1,319]:Integer, [-1,320]:Integer, [-1,321]:Integer, [-1,322]:Integer, [-1,323]:Integer, [-1,324]:Integer, [-1,325]:Integer, [-1,326]:Integer, [-1,327]:Integer, [-1,328]:Float@double, [-1,336]:Float@double, [-1,344]:Float@double, [-1,352]:Float@double, [-1,360]:Float@double, [-1,368]:Float@double, [-1,376]:Integer, [-1,377]:Integer, [-1,378]:Integer, [-1,379]:Integer, [-1,380]:Integer, [-1,381]:Integer, [-1,382]:Integer, [-1,383]:Integer, [-1,384]:Integer, [-1,385]:Integer, [-1,386]:Integer, [-1,387]:Integer, [-1,388]:Inte...

Farther down the error message says:

Cannot deduce type of insertvalue ins   %138 = insertvalue { { i64, i64, i64, i64, i64, i64, double, double, double, double, double, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, double, double, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, double, double, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] }, { { [2 x double], [2 x double], i64, i64 }, [1 x i64] } }, { {} addrspace(10)*, [3 x i64] }, {} addrspace(10)* } %134, { {} addrspace(10)*, [3 x i64] } %unbox2.i.unpack286435, 1, !dbg !20 size: 32 TT: {}

Caused by:
Stacktrace:
 [1] #60
   @ ./tuple.jl:461
 [2] afoldl
   @ ./operators.jl:545
 [3] filter_rec
   @ ./tuple.jl:461
 [4] filter
   @ ./tuple.jl:464
 [5] filter
   @ ./tuple.jl:0
within MethodInstance for filter(::Oceananigans.Fields.var"#87#91", ::Tuple{Field{Face, Center, Center, Nothing, RectilinearGrid{Float64, Periodic, Periodic, Bounded, Float64, Float64, Float64, OffsetArrays.OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}, OffsetArrays.OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}, OffsetArrays.OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}, CPU}, Tuple{Colon, Colon, Colon}, OffsetArrays.OffsetArray{Float64, 3, Array{Float64, 3}}, Float64, FieldBoundaryConditions{BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Flux, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Flux, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Flux, Nothing}}, Nothing, Oceananigans.Fields.FieldBoundaryBuffers{Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing}}, Field{Center, Face, Center, Nothing, RectilinearGrid{Float64, Periodic, Periodic, Bounded, Float64, Float64, Float64, OffsetArrays.OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}, OffsetArrays.OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}, OffsetArrays.OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}, CPU}, Tuple{Colon, Colon, Colon}, OffsetArrays.OffsetArray{Float64, 3, Array{Float64, 3}}, Float64, FieldBoundaryConditions{BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Flux, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Flux, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Flux, Nothing}}, Nothing, Oceananigans.Fields.FieldBoundaryBuffers{Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing}}, Field{Center, Center, Face, Nothing, RectilinearGrid{Float64, Periodic, Periodic, Bounded, Float64, Float64, Float64, OffsetArrays.OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}, OffsetArrays.OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}, OffsetArrays.OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}, CPU}, Tuple{Colon, Colon, UnitRange{Int64}}, OffsetArrays.OffsetArray{Float64, 3, Array{Float64, 3}}, Float64, FieldBoundaryConditions{BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Nothing, Nothing, BoundaryCondition{Oceananigans.BoundaryConditions.Flux, Nothing}}, Nothing, Oceananigans.Fields.FieldBoundaryBuffers{Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing}}})

And the stacktrace says:

Stacktrace:
  [1] #60
    @ ./tuple.jl:461 [inlined]
  [2] afoldl
    @ ./operators.jl:545 [inlined]
  [3] filter_rec
    @ ./tuple.jl:461 [inlined]
  [4] filter
    @ ./tuple.jl:464 [inlined]
  [5] filter
    @ ./tuple.jl:0 [inlined]
  [6] diffejulia_filter_21024_inner_1wrap
    @ ./tuple.jl:0
  [7] macro expansion
    @ ~/.julia/packages/Enzyme/BRtTP/src/compiler.jl:8137 [inlined]
  [8] enzyme_call(::Val{…}, ::Ptr{…}, ::Type{…}, ::Val{…}, ::Val{…}, ::Type{…}, ::Type{…}, ::Const{…}, ::Type{…}, ::Const{…}, ::MixedDuplicated{…}, ::Tuple{…}, ::Nothing)
    @ Enzyme.Compiler ~/.julia/packages/Enzyme/BRtTP/src/compiler.jl:7703
  [9] (::Enzyme.Compiler.AdjointThunk{Ptr{…}, Const{…}, MixedDuplicated{…}, Tuple{…}, 1, Nothing})(::Const{typeof(filter)}, ::Const{Oceananigans.Fields.var"#87#91"}, ::Vararg{Any})
    @ Enzyme.Compiler ~/.julia/packages/Enzyme/BRtTP/src/compiler.jl:7508
 [10] runtime_generic_rev(activity::Type{…}, runtimeActivity::Val{…}, width::Val{…}, ModifiedBetween::Val{…}, tape::Enzyme.Compiler.Tape{…}, f::typeof(filter), df::Nothing, primal_1::Oceananigans.Fields.var"#87#91", shadow_1_1::Nothing, primal_2::Tuple{…}, shadow_2_1::Base.RefValue{…})
    @ Enzyme.Compiler ~/.julia/packages/Enzyme/BRtTP/src/rules/jitrules.jl:682
 [11] #fill_halo_regions!#83
    @ ~/.julia/packages/Oceananigans/EYQtB/src/Fields/field_tuples.jl:69 [inlined]
 [12] diffejulia__fill_halo_regions__83_17571_inner_1wrap
    @ ~/.julia/packages/Oceananigans/EYQtB/src/Fields/field_tuples.jl:0
 [13] macro expansion
    @ ~/.julia/packages/Enzyme/BRtTP/src/compiler.jl:8137 [inlined]
 [14] enzyme_call
    @ ~/.julia/packages/Enzyme/BRtTP/src/compiler.jl:7703 [inlined]
 [15] AdjointThunk
    @ ~/.julia/packages/Enzyme/BRtTP/src/compiler.jl:7508 [inlined]
 [16] runtime_generic_rev(activity::Type{…}, runtimeActivity::Val{…}, width::Val{…}, ModifiedBetween::Val{…}, tape::Enzyme.Compiler.Tape{…}, f::Oceananigans.Fields.var"##fill_halo_regions!#83", df::Nothing, primal_1::@Kwargs{…}, shadow_1_1::Nothing, primal_2::typeof(Oceananigans.BoundaryConditions.fill_halo_regions!), shadow_2_1::Nothing, primal_3::@NamedTuple{…}, shadow_3_1::Base.RefValue{…}, primal_4::Clock{…}, shadow_4_1::Clock{…}, primal_5::@NamedTuple{…}, shadow_5_1::Base.RefValue{…})
    @ Enzyme.Compiler ~/.julia/packages/Enzyme/BRtTP/src/rules/jitrules.jl:689
 [17] fill_halo_regions!
    @ ~/.julia/packages/Oceananigans/EYQtB/src/Fields/field_tuples.jl:56 [inlined]
 [18] fill_halo_regions!
    @ ~/.julia/packages/Oceananigans/EYQtB/src/Fields/field_tuples.jl:0 [inlined]
 [19] diffejulia_fill_halo_regions__17431_inner_1wrap
    @ ~/.julia/packages/Oceananigans/EYQtB/src/Fields/field_tuples.jl:0
 [20] macro expansion
    @ ~/.julia/packages/Enzyme/BRtTP/src/compiler.jl:8137 [inlined]
 [21] enzyme_call
    @ ~/.julia/packages/Enzyme/BRtTP/src/compiler.jl:7703 [inlined]
 [22] AdjointThunk
    @ ~/.julia/packages/Enzyme/BRtTP/src/compiler.jl:7508 [inlined]
 [23] runtime_generic_rev(activity::Type{…}, runtimeActivity::Val{…}, width::Val{…}, ModifiedBetween::Val{…}, tape::Enzyme.Compiler.Tape{…}, f::typeof(Core.kwcall), df::Nothing, primal_1::@NamedTuple{…}, shadow_1_1::Nothing, primal_2::typeof(Oceananigans.BoundaryConditions.fill_halo_regions!), shadow_2_1::Nothing, primal_3::@NamedTuple{…}, shadow_3_1::Base.RefValue{…}, primal_4::Clock{…}, shadow_4_1::Clock{…}, primal_5::@NamedTuple{…}, shadow_5_1::Base.RefValue{…})
    @ Enzyme.Compiler ~/.julia/packages/Enzyme/BRtTP/src/rules/jitrules.jl:689
 [24] #update_state!#49
    @ ~/.julia/packages/Oceananigans/EYQtB/src/Models/HydrostaticFreeSurfaceModels/update_hydrostatic_free_surface_model_state.jl:41
 [25] update_state!
    @ ~/.julia/packages/Oceananigans/EYQtB/src/Models/HydrostaticFreeSurfaceModels/update_hydrostatic_free_surface_model_state.jl:32 [inlined]
 [26] update_state!
    @ ~/.julia/packages/Oceananigans/EYQtB/src/Models/HydrostaticFreeSurfaceModels/update_hydrostatic_free_surface_model_state.jl:29 [inlined]
 [27] #time_step!#8
    @ ~/.julia/packages/Oceananigans/EYQtB/src/TimeSteppers/quasi_adams_bashforth_2.jl:123
 [28] time_step!
    @ ~/.julia/packages/Oceananigans/EYQtB/src/TimeSteppers/quasi_adams_bashforth_2.jl:76 [inlined]
 [29] viscous_hydrostatic_turbulence
    @ ~/Projects/DifferentOcean/differentiated_hydrostatic_turbulence.jl:48 [inlined]
 [30] viscous_hydrostatic_turbulence
    @ ~/Projects/DifferentOcean/differentiated_hydrostatic_turbulence.jl:0 [inlined]
 [31] diffejulia_viscous_hydrostatic_turbulence_54612_inner_1wrap
    @ ~/Projects/DifferentOcean/differentiated_hydrostatic_turbulence.jl:0
 [32] macro expansion
    @ ~/.julia/packages/Enzyme/BRtTP/src/compiler.jl:8137 [inlined]
 [33] enzyme_call
    @ ~/.julia/packages/Enzyme/BRtTP/src/compiler.jl:7703 [inlined]
 [34] CombinedAdjointThunk
    @ ~/.julia/packages/Enzyme/BRtTP/src/compiler.jl:7476 [inlined]
 [35] autodiff
    @ ~/.julia/packages/Enzyme/BRtTP/src/Enzyme.jl:491 [inlined]
 [36] autodiff
    @ ~/.julia/packages/Enzyme/BRtTP/src/Enzyme.jl:537 [inlined]
 [37] autodiff(::ReverseMode{…}, ::typeof(viscous_hydrostatic_turbulence), ::Active{…}, ::Duplicated{…}, ::Const{…}, ::Const{…}, ::Const{…}, ::Const{…}, ::Const{…})
    @ Enzyme ~/.julia/packages/Enzyme/BRtTP/src/Enzyme.jl:504
 [38] top-level scope
    @ ~/Projects/DifferentOcean/differentiated_hydrostatic_turbulence.jl:78
 [39] include(fname::String)
    @ Base.MainInclude ./client.jl:489
 [40] top-level scope
    @ REPL[12]:1
in expression starting at /Users/gregorywagner/Projects/DifferentOcean/differentiated_hydrostatic_turbulence.jl:78
Some type information was truncated. Use `show(err)` to see complete types.

@glwagner
Copy link
Member Author

It's coming from this annoying function:

function fill_halo_regions!(maybe_nested_tuple::Union{NamedTuple, Tuple}, args...; kwargs...)
flattened = flattened_unique_values(maybe_nested_tuple)
# Sort fields into ReducedField and Field with non-nothing boundary conditions
fields_with_bcs = filter(f -> !isnothing(boundary_conditions(f)), flattened)
reduced_fields = filter(f -> f isa ReducedField, fields_with_bcs)
for field in reduced_fields
fill_halo_regions!(field, args...; kwargs...)
end
# MultiRegion fields are considered windowed_fields (indices isa MultiRegionObject))
windowed_fields = filter(f -> !(f isa FullField), fields_with_bcs)
ordinary_fields = filter(f -> (f isa FullField) && !(f isa ReducedField), fields_with_bcs)
# Fill halo regions for reduced and windowed fields
for field in windowed_fields
fill_halo_regions!(field, args...; kwargs...)
end
# Fill the rest
if !isempty(ordinary_fields)
grid = first(ordinary_fields).grid
tupled_fill_halo_regions!(ordinary_fields, grid, args...; kwargs...)
end
return nothing
end

which we can certainly simplify. We can also play around this avoiding hitting this line altogether, at least for a certain class of models.

@glwagner
Copy link
Member Author

Ok problem may be fixed.

@wsmoses
Copy link
Collaborator

wsmoses commented Oct 24, 2024

oh cool, well in any case this is why we should make PRs/land everything we want to experiment with [and things that fail I can fix]

@glwagner
Copy link
Member Author

oh cool, well in any case this is why we should make PRs/land everything we want to experiment with [and things that fail I can fix]

100% --- it seems like the problems here are all on the Oceananigans side

@jlk9
Copy link
Collaborator

jlk9 commented Oct 25, 2024

Is this correct in lines 274=275 of enzyme_test?

     x = y = (0, 2π)
     z = 1

When I run it locally it returns this error:

Got exception outside of a @test
  ArgumentError: z length(1) must be 2.

I found a few possible errors locally and can push.

@glwagner
Copy link
Member Author

glwagner commented Oct 25, 2024

Can you link to the code? Do you mean these lines:

x = y = (0, 2π)
z = 1
ν₀ = 1e-2
grid = RectilinearGrid(arch, size=(Nx, Ny, 1); x, y, z, topology=(Periodic, Periodic, Bounded))

Yes that's a bug -- fixed now

PS to generate a link to code, click the "..." on the left side. You'll get an option to create a permalink. You can edit the link to refer to a range of lines after creating it as I did above (to see what I wrote, attempt to edit my post and you will see the hyperlink)

image

@glwagner
Copy link
Member Author

Here's another code snipped that I have been using to test this:

using Oceananigans
using Oceananigans.Utils: with_tracers
using Random
using Enzyme

Random.seed!(123)
arch = CPU()
Nx = Ny = 32
x = y = (0, 2π)
z = (0, 1)
g = 4^2
c = sqrt(g)

grid = RectilinearGrid(arch, size=(Nx, Ny, 1); x, y, z, topology=(Periodic, Periodic, Bounded))
closure = ScalarDiffusivity=1e-2)
momentum_advection = Centered(order=2)
free_surface = ExplicitFreeSurface(gravitational_acceleration=g)
model = HydrostaticFreeSurfaceModel(; grid, momentum_advection, free_surface, closure)

ϵ(x, y, z) = 2randn() - 1
set!(model, u=ϵ, v=ϵ)

u_init = deepcopy(model.velocities.u)
v_init = deepcopy(model.velocities.v)

Δx = minimum_xspacing(grid)
Δt = 0.01 * Δx / c
for n = 1:10
    time_step!(model, Δt; euler=true)
end

u_truth = deepcopy(model.velocities.u)
v_truth = deepcopy(model.velocities.v)

function set_viscosity!(model, viscosity)
    new_closure = ScalarDiffusivity=viscosity)
    names = ()
    new_closure = with_tracers(names, new_closure)
    model.closure = new_closure
    return nothing
end

function viscous_hydrostatic_turbulence(ν, model, u_init, v_init, Δt, u_truth, v_truth)
    # Initialize the model
    model.clock.iteration = 0
    model.clock.time = 0
    #model.clock.last_Δt = Inf
    set_viscosity!(model, ν)
    #set!(model, u=u_init, v=v_init, η=0)
    set!(model, u=u_init, v=v_init)
    fill!(parent(model.free_surface.η), 0)

    # Step it forward
    for n = 1:10
        time_step!(model, Δt; euler=true)
    end

    # Compute the sum square error
    u, v, w = model.velocities
    Nx, Ny, Nz = size(model.grid)
    err = 0.0
    for j = 1:Ny, i = 1:Nx
        err += @inbounds (u[i, j, 1] - u_truth[i, j, 1])^2 +
                         (v[i, j, 1] - v_truth[i, j, 1])^2
    end

    return err::Float64
end

# Use a manual finite difference to compute a gradient
Δν = 1e-6
ν1 = 1.1e-2
ν2 = ν1 + Δν
e1 = viscous_hydrostatic_turbulence(ν1, model, u_init, v_init, Δt, u_truth, v_truth)
e2 = viscous_hydrostatic_turbulence(ν2, model, u_init, v_init, Δt, u_truth, v_truth)
Δe = e2 - e1
ΔeΔν = (e2 - e1) / Δν

@info "Finite difference computed: $ΔeΔν"

@info "Now with autodiff..."
start_time = time_ns()

# Use autodiff to compute a gradient
dmodel = Enzyme.make_zero(model)
dedν = autodiff(set_runtime_activity(Enzyme.Reverse),
                viscous_hydrostatic_turbulence,
                Active(ν1),
                Duplicated(model, dmodel),
                Const(u_init),
                Const(v_init),
                Const(Δt),
                Const(u_truth),
                Const(v_truth))

@info "Automatically computed: $dedν."
@info "Elapsed time: " * prettytime(1e-9 * (time_ns() - start_time))

@jlk9
Copy link
Collaborator

jlk9 commented Oct 25, 2024

Cool! Two other possible bugs I found:

rel_error = abs(dedν[1][3] - ΔeΔν) / abs(ΔeΔν)

I think should be dedν[1][1] instead. And

for j = 1:Ny, i = 1:Nx
err += @inbounds (u[i, j, 3] - u_truth[i, j, 3])^2 +
(v[i, j, 3] - v_truth[i, j, 3])^2
end

produces an out of bounds error with the z axis. I replaced 3 with end.

Happy to push those changes.

@glwagner
Copy link
Member Author

Cool! Two other possible bugs I found:

rel_error = abs(dedν[1][3] - ΔeΔν) / abs(ΔeΔν)

I think should be dedν[1][1] instead. And

I noticed that. Why was it [1][3] before?

We can't use end so let me fix that stuff now

@jlk9
Copy link
Collaborator

jlk9 commented Oct 25, 2024

Not exactly sure why it was [1][3], but if u_init was an active variable at one point then we we would want to get its derivative which is there.

@glwagner
Copy link
Member Author

Not exactly sure why it was [1][3], but if u_init was an active variable at one point then we we would want to get its derivative which is there.

Hm ok got it, it depends on the status of what's returned

@jlk9
Copy link
Collaborator

jlk9 commented Nov 7, 2024

Ok, so autodiff on the current commit at ν_0 = 1e-2 produces 0.0 as an output, which seems a little odd. But the FD comparison supports it. Here's the FD results across a range of step sizes:

(Δν, ΔeΔν) = (1.0e-5, -1.7467288660086175e-7)
(Δν, ΔeΔν) = (1.0e-6, -1.7467602081020574e-9)
(Δν, ΔeΔν) = (1.0e-7, -1.749279074332283e-11)
(Δν, ΔeΔν) = (1.0e-8, -1.9038931512874885e-13)
(Δν, ΔeΔν) = (1.0e-9, 1.0210771799161064e-14)
(Δν, ΔeΔν) = (1.0e-10, -7.313147905552346e-14)

However, since |0 - FD result| / |FD result| = 1.0, the test technically fails.

Before you were applying AD to ν1 = ν_0 + delta ν = 1e-2 + 1e-6. There the AD result is 0.004914, and central difference FD is 0.00491540, so the test passes.

@glwagner
Copy link
Member Author

glwagner commented Nov 7, 2024

Ok, so autodiff on the current commit at ν_0 = 1e-2 produces 0.0 as an output, which seems a little odd.

I believe this is because ν_0 is the global minimum, so the derivative is 0 there. It might be good to check a slightly different viscosity to get a non-trivial derivative

@glwagner
Copy link
Member Author

glwagner commented Nov 7, 2024

The CPU tests should be back up, so when you make another commit all the tests will run hopefully. It looks like the regression tests are passing so this should be ready to merge soon.

@glwagner glwagner requested a review from jlk9 November 7, 2024 17:10
@glwagner
Copy link
Member Author

glwagner commented Nov 7, 2024

@jlk9 you will have to give it an approval since I opened the PR

@jlk9
Copy link
Collaborator

jlk9 commented Nov 7, 2024

Ok, so autodiff on the current commit at ν_0 = 1e-2 produces 0.0 as an output, which seems a little odd.

I believe this is because ν_0 is the global minimum, so the derivative is 0 there. It might be good to check a slightly different viscosity to get a non-trivial derivative

Perfect, 0 derivative is expected for ν_0 then. I just tweaked the AD test so it's on ν1 = ν_0 + delta ν. Hopefully all tests will pass.

jlk9
jlk9 previously approved these changes Nov 7, 2024
@jlk9
Copy link
Collaborator

jlk9 commented Nov 7, 2024

Some of the distributed tests are failing, and I wonder if the changes to fill_halo_regions! are causing it. Looking at the stack traces for distributed cpu unit tests:

MethodError: no method matching fill_send_buffers!(::Tuple{OffsetArray{Float64, 3, Array{Float64, 3}}, OffsetArray{Float64, 3, Array{Float64, 3}}, OffsetArray{Float64, 3, Array{Float64, 3}}}, ::Clock{Float64, Float64}, ::RectilinearGrid{Float64, FullyConnected, Periodic, Periodic, Float64, Float64, Float64, OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}, OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}, OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}, Distributed{CPU, false, Partition{Int64, Nothing, Nothing}, Tuple{Int64, Int64, Int64}, Int64, Tuple{Int64, Int64, Int64}, Oceananigans.DistributedComputations.RankConnectivity{Int64, Int64, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing}, MPI.Comm, Vector{MPI.Request}, Base.RefValue{Int64}}}, ::Val{:bottom_and_top})
  The function `fill_send_buffers!` exists, but no method is defined for this combination of argument types.
  Closest candidates are:
    fill_send_buffers!(::OffsetArray, ::Oceananigans.Fields.FieldBoundaryBuffers, ::Any, ::Val{:bottom_and_top})
     @ Oceananigans /central/scratch/esm/slurm-buildkite/oceananigans-distributed/4309/oceananigans-distributed/src/Fields/field_boundary_buffers.jl:181
    fill_send_buffers!(::OffsetArray, ::Oceananigans.Fields.FieldBoundaryBuffers, ::Any, ::Val{:bottom})
     @ Oceananigans /central/scratch/esm/slurm-buildkite/oceananigans-distributed/4309/oceananigans-distributed/src/Fields/field_boundary_buffers.jl:158
    fill_send_buffers!(::OffsetArray, ::Oceananigans.Fields.FieldBoundaryBuffers, ::Any, ::Val{:north})
     @ Oceananigans /central/scratch/esm/slurm-buildkite/oceananigans-distributed/4309/oceananigans-distributed/src/Fields/field_boundary_buffers.jl:156

The stack trace goes through fill_halo_regions! in BoundaryConditions/fill_halo_regions.jl, then the method of fill_halo_event! inside of DistributedComputations/halo_communication.jl which is where fill_send_buffers is called. Should it go through fill_halo_event! in BoundaryConditions/fill_halo_regions.jl instead?

Looking at the function handles:

function fill_halo_event!(c, fill_halos!, bcs, indices, loc, arch, grid, args...; kwargs...)

function fill_halo_event!(c, fill_halos!, bcs, indices, loc, arch, grid::DistributedGrid, buffers, args...; async = false, only_local_halos = false, kwargs...)

The main difference in the handles is the buffer argument (and the grid being specified as distributed). Back at the stack trace:

... buffers::Clock{Float64, Float64}, args::@NamedTuple{ ...

we see that the Clock object is getting treated as a buffer instead of part of args.

The PR doesn't change fill_halo_regions.jl where this version of fill_halo_regions! is called. But it does change fill_halo_regions! in src/Fields/field_tuples.jl which is also part of the stack trace. Maybe the rewrite of fill_halo_regions! there affected how the clock argument is passed?

@jlk9 jlk9 dismissed their stale review November 7, 2024 22:20

Potential issues with distributed tests.

@glwagner
Copy link
Member Author

glwagner commented Nov 7, 2024

Some of the distributed tests are failing, and I wonder if the changes to fill_halo_regions! are causing it. Looking at the stack traces for distributed cpu unit tests:

MethodError: no method matching fill_send_buffers!(::Tuple{OffsetArray{Float64, 3, Array{Float64, 3}}, OffsetArray{Float64, 3, Array{Float64, 3}}, OffsetArray{Float64, 3, Array{Float64, 3}}}, ::Clock{Float64, Float64}, ::RectilinearGrid{Float64, FullyConnected, Periodic, Periodic, Float64, Float64, Float64, OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}, OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}, OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}, Distributed{CPU, false, Partition{Int64, Nothing, Nothing}, Tuple{Int64, Int64, Int64}, Int64, Tuple{Int64, Int64, Int64}, Oceananigans.DistributedComputations.RankConnectivity{Int64, Int64, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing}, MPI.Comm, Vector{MPI.Request}, Base.RefValue{Int64}}}, ::Val{:bottom_and_top})
  The function `fill_send_buffers!` exists, but no method is defined for this combination of argument types.
  Closest candidates are:
    fill_send_buffers!(::OffsetArray, ::Oceananigans.Fields.FieldBoundaryBuffers, ::Any, ::Val{:bottom_and_top})
     @ Oceananigans /central/scratch/esm/slurm-buildkite/oceananigans-distributed/4309/oceananigans-distributed/src/Fields/field_boundary_buffers.jl:181
    fill_send_buffers!(::OffsetArray, ::Oceananigans.Fields.FieldBoundaryBuffers, ::Any, ::Val{:bottom})
     @ Oceananigans /central/scratch/esm/slurm-buildkite/oceananigans-distributed/4309/oceananigans-distributed/src/Fields/field_boundary_buffers.jl:158
    fill_send_buffers!(::OffsetArray, ::Oceananigans.Fields.FieldBoundaryBuffers, ::Any, ::Val{:north})
     @ Oceananigans /central/scratch/esm/slurm-buildkite/oceananigans-distributed/4309/oceananigans-distributed/src/Fields/field_boundary_buffers.jl:156

The stack trace goes through fill_halo_regions! in BoundaryConditions/fill_halo_regions.jl, then the method of fill_halo_event! inside of DistributedComputations/halo_communication.jl which is where fill_send_buffers is called. Should it go through fill_halo_event! in BoundaryConditions/fill_halo_regions.jl instead?

Looking at the function handles:

function fill_halo_event!(c, fill_halos!, bcs, indices, loc, arch, grid, args...; kwargs...)

function fill_halo_event!(c, fill_halos!, bcs, indices, loc, arch, grid::DistributedGrid, buffers, args...; async = false, only_local_halos = false, kwargs...)

The main difference in the handles is the buffer argument (and the grid being specified as distributed). Back at the stack trace:

... buffers::Clock{Float64, Float64}, args::@NamedTuple{ ...

we see that the Clock object is getting treated as a buffer instead of part of args.

The PR doesn't change fill_halo_regions.jl where this version of fill_halo_regions! is called. But it does change fill_halo_regions! in src/Fields/field_tuples.jl which is also part of the stack trace. Maybe the rewrite of fill_halo_regions! there affected how the clock argument is passed?

i'll look into this

@simone-silvestri
Copy link
Collaborator

Some of the distributed tests are failing, and I wonder if the changes to fill_halo_regions! are causing it. Looking at the stack traces for distributed cpu unit tests:

I see you have changed tupled_fill_halo_regions!. Tupled halo filling is not supported in distributed computations, so in the source code we extend it here

function tupled_fill_halo_regions!(full_fields, grid::DistributedGrid, args...; kwargs...)
for field in full_fields
fill_halo_regions!(field, args...; kwargs...)
end
end

In this PR it seems like this extension has been bypassed, and the tests are trying to send halos of field tuples. Maybe it's a useful information to debug

@jlk9
Copy link
Collaborator

jlk9 commented Nov 10, 2024

function fill_halo_regions!(maybe_nested_tuple::Union{NamedTuple, Tuple}, args...; kwargs...)
flattened = flattened_unique_values(maybe_nested_tuple)
# Check to find grid:
for f in flattened
if !isnothing(boundary_conditions(f))
if !(f isa ReducedField) && (f isa FullField)
grid = f.grid
break
end
end
end
return tupled_fill_halo_regions!(flattened, grid, args...; kwargs...)
end
function tupled_fill_halo_regions!(fields, grid, args...; kwargs...)

Probably not the way to permanently fix it, but I've just added grid back to the handle of tupled_fill_halo_regions! so it matches the function's handles elsewhere. If the distributed tests now pass then this is it.

src/Fields/field_tuples.jl Outdated Show resolved Hide resolved
src/Fields/field_tuples.jl Outdated Show resolved Hide resolved
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

4 participants