Derivatives in initialization equations

[baggepinnen]

MTK requires that all derivatives that are referred to in initialization equations and the varmap are present also in the regular equations. In modelica, one can specify der(y) = 0 in initialization equations despite der(y) not appearing in any regular equation, which is nice because it's common to want to initialize at zero derivative of outputs of block components.

Here's a link to a simple modelica model that demonstrates this

https://playground.modelica.university/?model=ICBwYXJhbWV0ZXIgUmVhbCBlPTAuOCAiQ29lZmZpY2llbnQgb2YgcmVzdGl0dXRpb24iOwoKICBjb25zdGFudCBSZWFsIGVwcz0xZS0zICJTbWFsbCBoZWlnaHQiOwogIEJvb2xlYW4gZG9uZSAiRmxhZyB3aGVuIHRvIHR1cm4gb2ZmIGdyYXZpdHkiOwogIFJlYWwgaCAiSGVpZ2h0IjsKICBSZWFsIHYgIlZlbG9jaXR5IjsKICBSZWFsIHk7CmluaXRpYWwgZXF1YXRpb24KICBoID0gNy4wICJJbml0aWFsIGhlaWdodCI7CiAgZGVyKHkpID0gMC4wICJJbml0aWFsIGJ1dCBleHByZXNzZWQgaW4gdGVybXMgb2YgZGVyIjsKICAvL3YgPSAwLjAgIkluaXRpYWwgYnV0IGV4cHJlc3NlZCBpbiB0ZXJtcyBvZiBkZXIiOwogIGRvbmUgPSBmYWxzZTsKZXF1YXRpb24KICB2ID0gZGVyKGgpOwogIHkgPSBoOwogIGRlcih2KSA9IGlmIGRvbmUgdGhlbiAwIGVsc2UgLTkuODE7CiAgd2hlbiB7aDwwLGg8LWVwc30gdGhlbgogICAgZG9uZSA9IGg8LWVwczsKICAgIHJlaW5pdCh2LCAtZSooaWYgaDwtZXBzIHRoZW4gMCBlbHNlIHByZSh2KSkpOwogIGVuZCB3aGVuOwo%253D&report=VGhlIGNsYXNzaWMgTW9kZWxpY2EgImJvdW5jaW5nIGJhbGwiLiAgUHJlc3MgdGhlICJQbGF5IiBidXR0b24gdG8gYW5pbWF0ZS4KCkNoYW5nZSB0aGUgY29lZmZpY2llbnQgb2YgcmVzaXR1dGlvbiB0byBzZWUgaG93IGl0IGNoYW5nZXMgdGhlIGJvdW5jaW5nIG9mIHRoZSBiYWxsLgoKPHN2ZyBoZWlnaHQ9IjQwMCIgd2lkdGg9IjEwMCUiIHZpZXdib3g9IjAgMCA0MDAgNDAwIj4KICA8Y2lyY2xlIGN4PSJ7eyAxMDArbm93LnRpbWUqMjUgfX0iIGN5PSJ7eyAxMDArNyoyNSAtMjUqbm93LmggfX0iIHI9IjQwIiBzdHJva2U9ImJsYWNrIiBzdHJva2Utd2lkdGg9IjMiIGZpbGw9InJlZCIgLz4KICA8cmVjdCB4PSIwIiB5PSJ7ezE0MCs3KjI1fX0iIHdpZHRoPSI0MDAiIGhlaWdodD0iMTAiIHN0cm9rZT0iYmxhY2siLz4KPC9zdmc%252BCgo8Y2hhcnQgc2lnbmFscz0iaCI%252BPC9jaGFydD4%253D

The model is reproduced below in case the link stops working

  parameter Real e=0.8 "Coefficient of restitution";

  constant Real eps=1e-3 "Small height";
  Boolean done "Flag when to turn off gravity";
  Real h "Height";
  Real v "Velocity";
  Real y;
initial equation
  h = 7.0 "Initial height";
  der(y) = 0.0 "Initial vel but expressed in terms of der of output";
  //v = 0.0 "Initial vel";
  done = false;
equation
  v = der(h);
  y = h;
  der(v) = if done then 0 else -9.81;
  when {h<0,h<-eps} then
    done = h<-eps;
    reinit(v, -e*(if h<-eps then 0 else pre(v)));
  end when;

[AayushSabharwal]

structural_simplify should be able to recognize that D(y) is present in the system and handle it appropriately. To avoid having to specify constraints that are binding as initialization_eqs are, we can also search defaults. The user can remove the initial condition without having to re-simplify the system. For a nicer API, structural_simplify could take a higher_order_derivatives keyword which takes an array of higher order derivatives to compute.

[AayushSabharwal]

I'm not fully sure about this, but maybe we can handle this outside of structural simplification too? If y is an observed variable y ~ f(..) then D(y) ~ D(f(..)) and we should be able to use expand_derivatives to reduce this. Recursively applying this procedure to the RHS and combining with dummy_sub for known derivatives should just work?