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

Parameter splitting breaks runtime generated jacobians #3294

Open
david-hofmann opened this issue Jan 8, 2025 · 0 comments
Open

Parameter splitting breaks runtime generated jacobians #3294

david-hofmann opened this issue Jan 8, 2025 · 0 comments
Labels
bug Something isn't working

Comments

@david-hofmann
Copy link

🐞: Setting split=true breaks the call of a runtime generated Jacobian:

Minimal Reproducible Example 👇

using ModelingToolkit
using ModelingToolkit: t_nounits as t, D_nounits as D

sys = let
           @parameters σ ρ β a::Bool
           @variables x(t) y(t) z(t)

           eqs = [D(x) ~ a * σ * (y - x),
                  D(y) ~ x *- z) - y,
                  D(z) ~ x * y - β * z]

           @named sys = ODESystem(eqs, t)
           structural_simplify(sys; split=true)
end;

jac = generate_jacobian(sys, expression=Val{false})[1]
jac(unknowns(sys), parameters(sys), t)

Error & Stacktrace ⚠️

ERROR: BoundsError: attempt to access Tuple{Vector{SymbolicUtils.BasicSymbolic{Real}}, Vector{Any}, Num} at index [4]
Stacktrace:
 [1] getindex
   @ ./tuple.jl:31 [inlined]
 [2] macro expansion
   @ ~/.julia/packages/RuntimeGeneratedFunctions/M9ZX8/src/RuntimeGeneratedFunctions.jl:162 [inlined]
 [3] macro expansion
   @ ./none:0 [inlined]
 [4] generated_callfunc
   @ ./none:0 [inlined]
 [5] (::RuntimeGeneratedFunctions.RuntimeGeneratedFunction{…})(::Vector{…}, ::Vector{…}, ::Num)
   @ RuntimeGeneratedFunctions ~/.julia/packages/RuntimeGeneratedFunctions/M9ZX8/src/RuntimeGeneratedFunctions.jl:150
 [6] top-level scope
   @ REPL[42]:1
Some type information was truncated. Use `show(err)` to see complete types.

Maybe relevant to note that with split=false the generated Jacobian has 3 parameters:

RuntimeGeneratedFunction(#=in Symbolics=#, #=using Symbolics=#, :((ˍ₋arg1, ˍ₋arg2, t)->#= /home/david/.julia/packages/Symbolics/PxO3a/src/build_function.jl:340 =# @inbounds(begin
              #= /home/david/.julia/packages/Symbolics/PxO3a/src/build_function.jl:340 =#
              begin
                  #= /home/david/.julia/packages/SymbolicUtils/jf8aQ/src/code.jl:385 =#
                  #= /home/david/.julia/packages/SymbolicUtils/jf8aQ/src/code.jl:386 =#
                  #= /home/david/.julia/packages/SymbolicUtils/jf8aQ/src/code.jl:387 =#
                  begin
                      begin
                          begin
                              begin
                                  begin
                                      #= /home/david/.julia/packages/SymbolicUtils/jf8aQ/src/code.jl:480 =#
                                      (SymbolicUtils.Code.create_array)(typeof(ˍ₋arg1), nothing, Val{2}(), Val{(3, 3)}(), (*)((*)(-1, ˍ₋arg2[1]), ˍ₋arg2[2]), (+)((*)(-1, ˍ₋arg1[3]), ˍ₋arg2[3]), ˍ₋arg1[2], (*)(ˍ₋arg2[1], ˍ₋arg2[2]), -1, ˍ₋arg1[1], 0, (*)(-1, ˍ₋arg1[1]), (*)(-1, ˍ₋arg2[4]))
                                  end
                              end
                          end
                      end
                  end
              end
          end)))

whereas with split=true it has 4:

RuntimeGeneratedFunction(#=in Symbolics=#, #=using Symbolics=#, :((ˍ₋arg1, ˍ₋arg2, ˍ₋arg3, t)->#= /home/david/.julia/packages/Symbolics/PxO3a/src/build_function.jl:340 =# @inbounds(begin
              #= /home/david/.julia/packages/Symbolics/PxO3a/src/build_function.jl:340 =#
              begin
                  #= /home/david/.julia/packages/SymbolicUtils/jf8aQ/src/code.jl:385 =#
                  #= /home/david/.julia/packages/SymbolicUtils/jf8aQ/src/code.jl:386 =#
                  #= /home/david/.julia/packages/SymbolicUtils/jf8aQ/src/code.jl:387 =#
                  begin
                      begin
                          begin
                              begin
                                  begin
                                      #= /home/david/.julia/packages/SymbolicUtils/jf8aQ/src/code.jl:480 =#
                                      (SymbolicUtils.Code.create_array)(typeof(ˍ₋arg1), nothing, Val{2}(), Val{(3, 3)}(), (*)((*)(-1, ˍ₋arg3[1]), ˍ₋arg2[3]), (+)((*)(-1, ˍ₋arg1[3]), ˍ₋arg2[1]), ˍ₋arg1[2], (*)(ˍ₋arg3[1], ˍ₋arg2[3]), -1, ˍ₋arg1[1], 0, (*)(-1, ˍ₋arg1[1]), (*)(-1, ˍ₋arg2[2]))
                                  end
                              end
                          end
                      end
                  end
              end
          end)))

It would be nice if the Jacobian treated the parameters the same way as parameters(sys) does. If that is not possible then having a function that returns the parameters in their split structure would help.

Thank you!!

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
bug Something isn't working
Projects
None yet
Development

No branches or pull requests

1 participant