Minimum working example for ModelingToolkit's CTarget

I’m trying to use ModelingToolkit.build_function to generate a C function with “vector-like” inputs, and “matrix-like” outputs. I would prefer not to have the left-hand-side of the equations in the state vector, but I can’t figure out how to build a CTarget function without including the left-hand-side in the “variables” argument.

My Attempt

import Pkg
using ModelingToolkit

@variables x[1:2,1:2] y[1:4]
eqs = x[:] .~ y

build_function(eqs, [y..., x...]; 


"void diffeqf(double* du, double* RHS1) {\n  du[4] = RHS1[0];\n  du[5] = RHS1[2];\n  du[6] = RHS1[1];\n  du[7] = RHS1[3];\n}\n"

The output above has du first indexed at 4, because it is the 4th through 8th element in the state vector. I’d prefer for du and RHS1 to be entirely separate. How can I best accomplish this?

If I leave out the x... in the state vector argument to build_function

julia> @variables x[1:2,1:2] y[1:4]
julia> eqs = x[:] .~ y

julia> build_function(eqs, y, target=ModelingToolkit.CTarget())
ERROR: MethodError: no method matching +(::Nothing, ::Int64)
Closest candidates are:
  +(::Any, ::Any, ::Any, ::Any...) at operators.jl:538
  +(::ChainRulesCore.One, ::Any) at /home/joe/.julia/packages/ChainRulesCore/cpHLu/src/differential_arithmetic.jl:94
  +(::ChainRulesCore.Zero, ::Any) at /home/joe/.julia/packages/ChainRulesCore/cpHLu/src/differential_arithmetic.jl:63

It’s just not supported yet. CTarget would need more work.

1 Like

I might have time to look into this in 2021. I think all that would be required is some top-level parsing changes, since the equations produced are correct in this simple Matrix case, the indices just need to be offset.

Should I file an issue on GitHub? Or hold off for now, since this isn’t a supported case?

file an issue

Update: #736