Interpretation of x(time.start) for delay buffer initialization

[henrikt-ma]

What is the value of y in the simulation of this model?

model DelayInitialization
  Real x(start = 1.0, fixed = true);
  Real y = delay(x, 5.0);
  discrete Real z;
initial equation
  z = y;
equation
  der(x) = 0.0;
  when not initial() then
    reinit(x, 2.0);
  end when;
  when false then
    z = 0.0;
  end when;
  annotation(experiment(StopTime = 1.0));
end DelayInitialization;

In System Modeler, the value is 1.0, but it could be argued that the specification says it should be 2.0.

The reason I also included z in the model is to highlight that we at least need to make the interpretation that x(0.0) actually simplifies to x in the initialization problem, and for me it would seem natural that this value is also the value used to initialize the delay buffer.

What is the value of x in other tools?

[gkurzbach]

SimulationX: x=2, y=2, z=1

[henrikt-ma]

SimulationX: x=2, y=2, z=1

Oh no. This will make it harder to define our result as the correct one. :(

[perost]

SimulationX: x=2, y=2, z=1

Same in OpenModelica.

[MarkusOlssonModelon]

Modelon Impact: x=2, y=1, z=1

[HansOlsson]

Dymola (and 3D Experience platform): x=2, y=2, z=1

[henrikt-ma]

Summarizing the tool reports:

• y = 1: System Modeler, Modelon Impact
• y = 2: Simulation X, OpenModelica, Dymola

So there is a majority with y = 2, but I will try to explain why I find this strange.

When defining delay(x, 5.0) as x(0) for time - 5.0 ≤ time.start, I like to think of this as the delay memory for x being initialized with the constant value x(0) for the entire history up to time.start. So, with x' representing the delay memory, we have

x' = if time ≤ time.start then x(0) else x

and then delay(x, 5.0) = x'(time - 5.0) is just about reading past values from the delay memory.

In our example, x is 2 from time 0 and forward, so the delay memory can be simplified:

x' = if time ≤ time.start then x(0) else 2

Also, z = 1 (upon which all tools agree) shows that during initialization we can see that x'(-5.0) = 1, so the constant value in the delay buffer up to time.start must be 1:

x' = if time ≤ time.start then 1 else 2

With this interpretation, y is 1 until time 5, when it becomes 2.

If y would be 2 throughout the simulation, I can't see how one can rescue the mental picture of delay(x, 5.0) as a time-shifted variant of an extrapolated x. So while that is a drawback of going with y = 2, are there also some clear advantages?

[casella]

The MLS says:

evaluates to 𝑒𝑥𝑝𝑟(time - 𝑑𝑒𝑙𝑎𝑦𝑇𝑖𝑚𝑒) for 𝚝𝚒𝚖𝚎>𝚝𝚒𝚖𝚎.𝚜𝚝𝚊𝚛𝚝+𝑑𝑒𝑙𝑎𝑦𝑇𝑖𝑚𝑒 and 𝑒𝑥𝑝𝑟(time.start) for 𝚝𝚒𝚖𝚎≤𝚝𝚒𝚖𝚎.𝚜𝚝𝚊𝚛𝚝+𝑑𝑒𝑙𝑎𝑦𝑇𝑖𝑚𝑒.

So, the question is, what is the value of expr(time.start), since there are events happening at that time point? Is it before the events take place or after? When we plot results, usually we plot the value after the last event has been processed. But that's not necessarily the best idea. In fact, one could argue that for initialization one should try to "look into the past" as much as possible, so maybe a better interpretation could be to consider expr(time.start) before any event has taken place there.

No strong opinions. For sure this should be clarified, one way or the other.

[HansOlsson]

After some quick thinking (still catching up after the conference, but also putting out fires, etc) I think there's a good reason to prefer the parameter-value x.start.

Users normally expect delay to break loops (see also #3726 ), but if delay(x, ...) should equal the initial value of x there's a direct dependency during initialization. This will occur both if x depends on this delayed value, and also if der(x) depends on it and the initialization links x and der(x); in all of those cases the initialization may be more complicated and error-prone.

For the future we might add more possibilities for initializing the delay-buffer, see #3633, so to me it makes sense to use x.start for the time being.

[henrikt-ma]

I don't think the "parameter-value x.start" is a feasible solution, since the delayed expression can be more complicated than just a variable. Also, when using delay to break loops, I would generally not expect the variable u in delay(u, …) to be determined by a start-attribute, but to be determined by the equations of the component producing u.

However, if u is a variable with fixed initialization (u = u.start), there will be decoupling at u in the initialization problem, also when delay(u, …) is simplified to just u. Whether this can be meaningfully utilized without running into overdetermined initialization is another question…

For now, I think simplifying delay(u, …) do u during initialization is as close to u.start as we can get.

Then I agree that giving the user control over the initialization would be a great feature for breaking loops also in the initialization problem. It would seem natural to me to borrow ideas from spatialDistribution:

delay(u, 1.0, initialTimes = {-1.0, -0.5, -0.5, 0.0}, initialValues = {7.0, 9.0, 8.0, 8.5})

[HansOlsson]

I don't think the "parameter-value x.start" is a feasible solution, since the delayed expression can be more complicated than just a variable. Also, when using delay to break loops, I would generally not expect the variable u in delay(u, …) to be determined by a start-attribute, but to be determined by the equations of the component producing u.

However, if u is a variable with fixed initialization (u = u.start), there will be decoupling at u in the initialization problem, also when delay(u, …) is simplified to just u. Whether this can be meaningfully utilized without running into overdetermined initialization is another question…

The decoupling when using the start-value comes automatically, and doesn't introduce a separate initialization problem for the delay.

There has been considerable push from users for ensuring that delay breaks algebraic loops, so I don't see any benefit of standardizing something that un-breaks during initialization without a clear physical background.

For now, I think simplifying delay(u, …) do u during initialization is as close to u.start as we can get.

To me simplifying it to u.start is simpler in practice and avoids the issues.

I know that in general delay takes an expression, but in practice it is often just a variable (if it is an expression, just take the start-values of the variables). Easy, pragmatic, and can later be extended.

Then I agree that giving the user control over the initialization would be a great feature fore breaking loops also in the initialization problem. It would seem natural to be to borrow ideas from spatialDistribution:

delay(u, 1.0, initialTimes = {-1.0, -0.5, -0.5, 0.0}, initialValues = {7.0, 9.0, 8.0, 8.5})

Something like that, but as I see it the main issue is initialization when continuing simulations - which is more of a tool issue; but they go hand in hand.