弹跳球

对弹跳球的建模

前面的例子里,我们看到了某些事件和时间的关系。这些所谓“时间事件”(time event)只是其中一种“事件”(event)。我们将在本节中研究另一类事件,即状态事件(state event)。状态事件是依赖于解的轨迹的一类事件。

处理状态事件相对而言更为复杂。不同于能够先验地知道发生时间的时间事件,状态事件依赖于解的轨迹。因此,我们无法完全避免去“搜索”发生事件的时点。

为了看到状态事件的效果,我们考虑一个弹跳球在平坦的水平面上弹跳时的行为。球在平面以上时会在重力作用下加速。而球最终与平面接触时,它会遵从以下关系从表面反弹:

\[v_\text{final} = -e v_\text{initial}\]

其中\(v_\text{final}\)为球在碰撞瞬间(在垂直方向)的速度\(v_\text{initial}\)为球在碰撞前的速度。而\(e\)则是恢复系数。恢复系数是对球在碰撞后所保留动量占比大小的度量。

如果我们把以上内容转化为Modelica,可能会是以下的样子:

model BouncingBall "The 'classic' bouncing ball model"
  type Height=Real(unit="m");
  type Velocity=Real(unit="m/s");
  parameter Real e=0.8 "Coefficient of restitution";
  parameter Height h0=1.0 "Initial height";
  Height h;
  Velocity v;
initial equation
  h = h0;
equation
  v = der(h);
  der(v) = -9.81;
  when h<0 then
    reinit(v, -e*pre(v));
  end when;
end BouncingBall;

在这个例子中,我们采用参数h0指定球在平面上的初始高度,参数e指定恢复系数。变量hv分别表示高度以及垂直速度。

本例有趣的地方在于模型中的公式。特别是其中的when语句:

equation
  v = der(h);
  der(v) = -9.81;
  when h<0 then
    reinit(v, -e*pre(v));
  end when;

when(当)语句由两部分组成。第一部分是条件表达式,表示事件发生的时刻。在本例中,事件发生在“当”高度h第一次降低到0以下时。when表达式的第二部分描述了在事件出现时所发生的事情。在本例中,v的值被reinit(重新初始化)操作符重新初始化了。我们可以用reinit操作符为状态设定新的初始化条件。从概念上讲,你可以认为reinit当成是在仿真过程中插入的initial equation。但是,reinit只会改变一个变量的值,而且这种改变总是显式的(也就是说它没有initial equation那么灵活)。在本例中,reinit语句会将v的值重新初始化,使之与在碰撞前的v(也就是pre(v))方向相反,大小变为原来的e倍。

假设h0为正值,重力不懈的拉扯确保球最终将贴近平面。在h0为1.0的情况下运行仿真,我们可以看到模型有以下行为:

../../../_images/BB1.png

在此图中,我们看到在0.48秒左右,球第一次和表面发生了碰撞。这是因为条件h<0在该时刻首次为真。需要注意的是,是什么使之成为状态事件(这个条件表达式引用了其他变量而不是变量time,而非我们在前面冷却例子内见到的样例。)

因此,仿真进行时假定球是在自由下落。一直到它发现在解的轨迹中某个特定的时间步长内,条件表达式h<0的值发生变化。倘若出现了这样的时间步长,求解器必须确定条件式的值变真的确定时刻。一旦该时刻被确定,求解器通过处理在该时间内影响了系统状态的when语句(如任何reinit语句)以计算该系统的状态。然后求解器以计算出的状态为初值重启积分计算。本例弹跳球内reinit语句的作用是,计算一个碰撞后让球重新(开始)上升的新值v

但要记住的是,一般情况下大多数的Modelica模型的解是由数值方法得出的。正如我们将很快看到,这会对我们考虑的离散行为产生深远的影响。这是因为,在所有的事件(时间或状态事件)的核心都是条件表达式,例如本例中的h<0

倘若我们模拟弹跳球的时间稍长,这种含义就会变得很清楚。在这种情况下,大多数Modelica语言的工具将提供这样的一个解:

../../../_images/BB2.png

显然,我们在看到轨迹后马上知道有些问题了。但到底是什么有问题呢?

数值精度

正如我们前面的暗示,答案就在于对条件表达式h<0的数值处理上。由于数值精度的问题,我们无法得知我们在开始计算本步时是到底在一个刚刚发生的事件之后,还是在事件即将发生之时。

为了解决这个问题,我们必须引入一定量的迟滞(死区)。这意味着在这种情况下,一旦条件h<0为真,我们必须离状态“足够远”才能让事件再次发生。换句话说,当h小于零时会发生事件。但在我们再次触发事件前,我们要求h必须大于某个\(\epsilon\)。换句话说,h并不是仅仅简单地大于零,h必须大于\(\epsilon\)(其中\(\epsilon\)是由求解器通过检查一系列缩放因子得到的)。

之前仿真里的问题是,每次球反弹时,h的峰值都会有些下降。这里峰值是指,球即将开始再次下降时h的值。最终,h峰值不足以超过该临界值\(\epsilon\)。反过来,这意味着when语句从未触发以及reinit语句不会再重置v。结果是球在继续不停地自由下落

所以一个显然的问题就是,如何实现我们真正的目标行为(球永远不会降到地面以下)。因此,我们必须对模型进行如下的一些小改动,:

model StableBouncingBall
  "The 'classic' bouncing ball model with numerical tolerances"
  type Height=Real(unit="m");
  type Velocity=Real(unit="m/s");
  parameter Real e=0.8 "Coefficient of restitution";
  parameter Height h0=1.0 "Initial height";
  constant Height eps=1e-3 "Small height";
  Boolean done;
  Height h;
  Velocity v;
initial equation
  h = h0;
  done = false;
equation
  v = der(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;
end StableBouncingBall;

应该指出的是,许多方法都可以用以解决这个问题。这里提出的解决方案仅是其中之一。在这种方法中,我们实际上创建了两个表面。其中一个在0的高度上,另一个在一个高度-eps(略低于0)。球“正常地”弹起时只会触发when语句的首个条件。然而,如果球在触地后不能足够高地回弹,从而“穿过”中的第一表面。我们检测到这点(以及球已经穿过的事实),并设置done标志符。done标志符的效果实际上就是关闭重力。

注意这里when语句的语法:

  when {h<0,h<-eps} then
    done = h<-eps;
    reinit(v, -e*(if h<-eps then 0 else pre(v)));
  end when;

特别是要注意,这个语句不只有一个条件表达式,而是两个。更准确地说,它实际上有向量形式的条件表达式。在本书的后面我们会介绍向量与数组,但在本章我们必须指出when既可以包含一个标量条件表达式,也可以包含矢量条件表达式。

如果一个when语句包含矢量条件表达式,那么当任意一个条件表达式为真时便会触发when语句内的语句。请仔细注意这断解释的语法。人们很有可能把如下形式Modelica代码:

when {a>0, b>0} then
  ...
end when;

理解为“当a大于零b大于零”。非常重要的是要避免一个十分常见的错误,即误认为以下两个when语句是等价的:

when {a>0, b>0} then
  ...
end when;

when a>0 or b>0 then
  ...
end when;

这两者并不等价的。为了理解上述差异,我们把条件表达式改为如下形式:

when {time>1, time>2} then
  ...
end when;

when time>1 or time>2 then
  ...
end when;

还记得我们之前对when语句中矢量表达式的陈述:当任意一个条件表达式为真时便会触发when语句内的语句。假设我们从time= 0开始运行仿真,直到time=3,那么when语句:

when {time>1, time>2} then
  ...
end when;

会被两次触发。第一次在当time>1为真时。另一次在time>2为真时。与此相反,在这种情况下:

when time>1 or time>2 then
  ...
end when;

仅存在一个条件表达式,且这个表达式仅仅一次变为真(当time>1为真﹍﹍并保持为真)。or操作符实际上屏蔽了第二项条件time>2,使得第二项条件在此情况下甚至可以不存在。换句话说,该条件仅一次变为真。其结果是,在when语句中的语句仅被触发一次。

关键是要记住的是,对于when语句,向量条件表示任何,而非逻辑或。此外,该语句仅在条件变为真的瞬间被激活。我们会在本章的后面部分讨论对比if和when之间的重要区别时,进一步探讨最后这句话的含义。