Generalized affects in MTK - what's a good interface?

Specific Domains Modelling & Simulations
1

Hi,

I’ve recently started a pull request in an attempt to add more generalized affect! support to MTK (and also to support discrete_events, which is taken on by others)

The reason I think this is needed is that in some instances, there’s a need to interface with a complex component that can’t be modeled easily by a set of differential equations. For example, a PV module could be controlled by a complex MPPT.

Anyway, I’m stuck on the right interface.

Currently, one can define an event as follows:

ODESystem(eqs, t; continuous_events = [ x ~ 0 ] => [v ~ -v])

Now, for a generalized affect! we can do something like:

ODESystem(eqs, t; continuous_events = [ x ~ 0 ] => affect!)

where affect! is a function. However, note that in this case the user doesn’t pass the variables/parameters of interest.

Therefore, with such an interface, MTK will have to expose all variables/parameters belonging to the ODESystem (note that if the affect! is defined as part of a sub-system, in principle, it could access all states and parameters of the entire system. This would break modularity, so let’s not go there).

For example, we could have:

function affect!(t, sts, pars)
   sts.v = -sts.v
end

However, MTK aggressively prunes unnecessary equations and variables. It can deduce (although it currently doesn’t) what variables that appear in affect equations should be “protected” from pruning. However, for the opaque affect! above, this doesn’t seem possible.

Another issue is that of “context”: to be useful, the affect! function should have access to some context/state it can interact with (e.g. the MPPT component).

The user could use closure to achieve this:

function make_affect(context)
   let context = context
      function (t,u,p)
         update!(context, u)
         p.C = new_c(context)
      end
   end
end

affect! = make_affect(MyContext())

However, this looks to me like an hack (and also a lot to expect from a library user).

So, an alternative interface could be for the user to pass not just the function, e.g.:

ODESystem(eqs, t; continuous_events = [ x ~ 0 ] => (affect!, context, (v = v, r = resistor.r), [C]))

This would allow MTK to pass the requested variables/parameters via kwargs:

function affect!(t, ctx; v, r, C)
   v = -C*v
   update_context!(ctx, r)
end

The main drawback I see in this solution is it’s cumbersome, but it’s the best I got.

What do you think?

DD

2

This is probably a better discussion for the MTK issues. You won’t catch Yingbo on here often.

Though I wonder if it should just be

ODESystem(eqs, t; continuous_events = [ x ~ 0 ] => affect!(:integrator,x,y))

which generates a function:

affect!(integrator,x,y)

where x is the integer for the index of x. That way you do integrator.u[x].

But honestly, that’s then getting very close to just writing it yourself with varmap_to_vars. It might be simpler to just add symbolic indexing support to AbstractVectorOfArray and then integrator.u[x] would just work on the symbolic x, or just make integrator[x] work on the symbolic x. Then writing a the callback by hand is not much harder than this other interface.

1 Like
3

Encapsulating the affect construction as a function call could work.

However, I don’t like the idea of exposing integrator this way: this breaks modularity.
I don’t want code simulating a bolt having access to parameters of the jet it is part of.

I think that:

function affect!(;x,y)
   x = x*y
end

is better than:

function affect!(integ; x, y)
   integ[x] = integ[x]*integ[y]
end

I think this can be done at zero-cost.
An alternative (which I don’t like) would be:

function affect!(t, sts, pars)
   sts.x = sts.x * sts.y
end

This is actually not far from passing the integrator - or rather a view of the integrator with only the variables requested in the call to affect! exposed (which now I think about it, is probably what you meant).

DD

4

Hmm…

the syntax I’ve suggested won’t work:

   x = x*y

will replace x. I’m pretty sure I can’t overload =.
Of course, I could use a similar syntax to Symbolics:

   x ~ x*y

But this will only confuse the users.

Maybe sts.x = sts.x * sts.y or sts[:x] = sts[:x] * sts[:y] is not that bad…

5

What about

ODESystem(eqs, t; continuous_events = [ x ~ 0 ] => (affect!, (v = v, r = resistor.r), [C]))

with

function affect!(integrator; v, r, C)
    integrator.u[v] = -integrator.p[C] * integrator.u[v] - integrator.u[r]
end

?

But honestly, that’s then getting very close to just writing it yourself with varmap_to_vars. It might be simpler to just add symbolic indexing support to AbstractVectorOfArray and then integrator.u[x] would just work on the symbolic x, or just make integrator[x] work on the symbolic x. Then writing a the callback by hand is not much harder than this other interface.

Working with Symbol would put strain on the Julia compiler and also ambiguous. For instance, what if x is a symbolic array? If x is namespaced deep in a component, the user has to write \_+ as well. It also kind of defeats the point of a-causal modeling as one needs to rewrite affect! when the model structure changes.