从定时自动机到极大极小代数：计算系统吞吐量与瓶颈

1. 为什么要专门给DES加一个时间维度第四章的定位先说一个我自己的感受。离散事件系统这本教材的前三章学起来其实很舒服自动机、语言、Petri网模型直观、概念有画面感。但等我把一条两机床流水线的自动机模型画好兴致勃勃地想回答这条线一小时到底能出多少个零件、瓶颈到底在哪台机器时发现完全使不上劲。自动机模型回答得了某个事件序列能不能发生回答不了什么时候发生更回答不了每隔多久发生一次。第四章解决的就是这个时间缺口。它给DES引入时间维度核心可以压缩成两条线一条是模型层的给自动机配上计时器形成定时自动机给Petri网配上时延形成定时事件图另一条是分析层的把一堆谁等谁的约束写成极大极小代数max-plus algebra下的方程组然后像普通线性代数那样求矩阵幂、特征值、特征向量把周期时间和吞吐量算出来。我自己把这章定义为从定性走向定量的分水岭。前几章学的是系统可能呈现哪些行为这一章学的是系统以什么节拍运行、最快能多快、瓶颈在哪。后面很多内容——随机定时自动机、马尔可夫链、排队网、监督控制——都要先默认你已经理解了这里的时间建模方式。直接跳过去读后面大概率会出现每个词都认识、但不知道在算什么的窘境。这篇笔记按原书第四章的主线整理我把教材的推导、课堂理解和跑代码验证结合在一起尽量用大白话把几个关键的弯弯绕绕讲透想给正在啃这本书的人省点功夫。2. 定时自动机与定时事件序列状态里多了一个时钟2.1 普通自动机到底缺了什么普通确定型自动机通常用事件集合、状态集合、状态转移函数、活跃事件函数和初始状态来刻画。它的运转逻辑是处在某个状态时有一批事件是使能的触发其中一个就跳到下一个状态。整个过程里没有任何变量代表时间。比如一台只有空闲和加工两个状态的机床在空闲时开始加工使能进入加工状态后加工完成使能再回到空闲。这个模型把顺序约束表达得很干净但它描述不了加工耗时2小时这件事。对只关心事件顺序的验证任务这足够对产能、节拍、瓶颈这类问题它缺一个最基本的维度。2.2 定时自动机的核心机制事件寿命与状态时钟原书引入定时自动机的方式是在普通自动机的状态上叠加一个时钟。系统停留在一个状态x时每个使能事件e都带有一个寿命参数可以理解为在当前状态下事件e最早还要等多久才能发生。真正的系统状态变成(x, v)其中v是已经在这个状态里停留的时间。当v增长到某个使能事件的寿命时该事件触发状态跳转v清零重新开始计时。和我最初预想的不一样定时自动机不是简单地在转移弧上标一个时间而是建模在当前状态下的等待时间。这两者的区别很实际同一条转移不同状态下可以有不同等待时长。比如超时关闭这个事件在机床空闲状态下寿命可能是300秒如果机床正在加工这个事件根本不使能。把寿命挂在状态上才能表达这种依赖关系。2.3 定时事件序列与dater函数换个视角记事件系统运行起来我们会记录下每个事件的触发时刻得到形如(e1,t1),(e2,t2),...的定时事件序列。这个序列本身很有用但它是一长条流水账很难直接做数学处理。原书在这里做了一个关键视角转换按事件类型拆开。具体来说对每一种事件i定义dater函数x_i(k)为事件i第k次发生的时刻。于是原来的一条定时事件序列变成若干条某类事件第k次发生的时刻数列。这个转换的重要意义在于制造系统、交通系统里同一个事件第k次和第k1次之间的间隔往往只取决于几条明确的等待约束。用数列表达以后下一步才能进入代数领域。dater这个称呼在不同教材里略有差异有的也叫计时函数意思就是给每个事件编号并记录发生时刻。2.4 一个超时模型的例子举个具体例子。设一台空闲机床在等待零件到达同时有一个超时关机事件。零件到达事件在空闲状态下的寿命不定取决于上游什么时候送来超时关机事件的寿命固定为5分钟。系统进入空闲状态的瞬间时钟从0开始走如果零件先到触发到达事件如果时钟先走到5触发超时关机。这个模型表达的逻辑是5分钟内没来零件就自动关掉在普通自动机里根本写不出来因为普通自动机没有已经等待了多少时间的概念。这一节能提供的直观就是时间信息应该放在状态停留过程里而不是当成一个独立事件塞进事件集合。理解了这个后面把时间约束翻译成dater方程时就不容易跑偏。3. 极大极小代数max与加法如何构成一门代数3.1 两条规则加法变成取最大值乘法变成加法极大极小代数的全部公理可以浓缩成两个运算a⊕b max(a,b)a⊗b ab。第一次看到会觉得很没必要——这不就是把取最大和做加法包装了一下吗但你把它们放进代数结构里看就会明白它为什么能成立一门代数。它满足分配律max(a,b) c max(ac, bc)这正对应普通代数里乘法对加法分配。它有零元ε-∞因为max(a, ε)a有单位元0因为a0a。用一张表对比普通代数和极大极小代数会更清楚概念普通代数极大极小代数基本操作ab, a×ba⊕bmax(a,b), a⊗bab零元0ε -∞单位元10矩阵乘法先乘后加、按项求和先加后取max特征值系统增长的幅度每个周期的时延增量为什么偏偏是max和加法因为DES的时间问题本质上是等待与耗时一个事件要等多条前驱路径都完成取最晚的那个时刻这就是max一段处理过程消耗的时长这就是加法。普通代数擅长描述数量的累积极大极小代数擅长描述时间延迟的累积两者结构高度平行所以整套线性代数工具都能搬过来用。3.2 矩阵乘法就是最长路径极大极小代数下的矩阵乘法可以写成(A⊗B)_ij max_k (A_ik B_kj)。把它翻译成图语言味道完全不同把矩阵A看成带权图A_ij就是从节点j到节点i的边的权值那么(A⊗A)_ij的值就是从j出发走两步到达i的所有路径里权值总和最大的一条更一般地A^m的(i,j)项就是从j到i的长度为m的路径中总延迟最大的一条。这意味着当我们把系统写成 x(k) A⊗x(k-1) 这样的递推时A^k的每一项不是在算一条固定路径而是把所有可能路径里最慢的那条挑出来。物理上这完全合理系统的某个事件第k次发生取决于所有导致它的环节里最晚完成的那一个。取max不是数学上的拍脑袋而是时间约束的本质。3.3 为什么这里敢说线性我刚开始学时有个长期困惑max明明是个强非线性运算教科书怎么敢说这是线性系统后来想明白了这里的线性不是数值意义上的而是代数结构意义上的。只要一个系统能写成状态向量 系数矩阵 ⊗ 状态向量的形式并且矩阵运算满足结合律、分配律那么矩阵幂、特征值、特征向量、线性方程组的求解逻辑就全部成立。这样说有点抽象换个类比。普通代数里你存了一笔钱第二年本息 a×本金在极大极小代数里一条产线第k个零件的完工时刻等于上一批零件某些时刻乘上加工时长再按取max的方式加和。两者形式上一模一样只是运算符号换了。所以极大极小系统在它自己的代数世界里就是线性系统这正是它的威力所在调度和产能问题可以被当作线性代数来处理。4. 系统矩阵与特征值吞吐量是怎么被算出来的4.1 两台串联机床的dater方程与系统矩阵前面讲的都是概念这一节落一个完整的算例。设一条流水线有两台机床A和B串联中间缓冲区最多放1个零件。机床A加工一个零件耗时2个单位机床B耗时3个单位。每个零件必须先在A上加工再进缓冲区再上B加工最后离开。记a_k为第k个零件在A上开工的时刻b_k为第k个零件在B上开工的时刻c_k为第k个零件完工离开的时刻。约束有四条。A一次只能做一个零件a_k ≥ a_{k-1}2B一次只能做一个零件b_k ≥ b_{k-1}3零件必须先A后Bb_k ≥ a_k2缓冲区容量为1A要开工第k1个零件必须等B已经把第k个零件取走a_{k1} ≥ b_k也就是a_k ≥ b_{k-1}。取每条约束的等号得到最早调度方程a_k max(a_{k-1}2, b_{k-1}) b_k max(b_{k-1}3, a_k2) c_k b_k3为了把b_k也写成只依赖上一时刻状态的形式做一次代换a_k2 max(a_{k-1}4, b_{k-1}2)所以 b_k max(b_{k-1}3, a_{k-1}4, b_{k-1}2) max(b_{k-1}3, a_{k-1}4)。同理 c_k max(b_{k-1}6, a_{k-1}7)。于是可以写成 x(k) A⊗x(k-1)其中 x(k)[a_k, b_k, c_k]^T矩阵A各行的含义如下x(k)a_{k-1}b_{k-1}c_{k-1}a_k20εb_k43εc_k76ε这里ε-∞表示对应位置没有依赖关系。我写过一段简单代码验证这个递推A_time 2 B_time 3 N 8 a [float(-inf)] * (N1) b [float(-inf)] * (N1) c [float(-inf)] * (N1) for k in range(1, N1): if k 1: a[k] 0 else: a[k] max(a[k-1] A_time, b[k-1]) b[k] max(b[k-1] B_time, a[k] A_time) c[k] b[k] B_time for k in range(1, N1): print(k, a[k], b[k], c[k])运行结果第1个零件 a0, b2, c5第2个 a2, b5, c8第3个 a5, b8, c11之后每3个时间单位出一件。周期时间λ3吞吐量就是1/3件/单位时间。4.2 特征值与关键环瓶颈藏在图里现在来看矩阵A对应的带权图。节点有三个a、b、c。边有a→a权2b→a权0a→b权4b→b权3a→c权7b→c权6。c没有出边所以它不参与任何环。需要找的是所有环里平均权值最大的一条a→a平均2b→b平均3a→b→a平均(40)/22。最大平均环权是3这就是极大极小特征值λ3。特征值在这里的物理含义是每推进一个零件周期时间前进多少也就是系统的周期时间。再看特征向量。解 A⊗v λ⊗v也就是让每个分量都加上λ3之后等式成立可取v[0,3,6]^T。三个分量的相对偏差3和6对应稳态下a、b、c三类事件的错位B总是在A开工后3个时间单位才开始C又在B开工后3个时间单位才完工。这些偏移不是随便定的它们正好让缓冲区约束和单机约束同时饱和。我一开始觉得特征向量这东西在实际调度里没什么用后来才发现它对应的是各事件的稳定相位关系。即使初始状态差得很远系统经过一段时间的追赶后三类事件的间隔也会收敛到这组相对偏移上来。理解这个排产时就知道该把哪类事件提前、哪类事件压后。4.3 改进哪台机器才有用临界环的决策价值有了特征值就能做一个很实际的决策实验如果把机床A的加工时间从2改成1会怎样我把加工参数改一下重算周期时间依然是3吞吐量纹丝不动。原因在图里看得很清楚决定λ3的临界环是b→b这条自环它只依赖B的加工时间3。A再快B每3个单位才能完成一个零件产量就是上不去。反过来如果把B的时间从3降到2周期时间立刻变成2吞吐量变成1/2。计算场景对比如下场景A加工时间B加工时间周期时间λ吞吐量原方案2331/3只加快A1331/3只加快B2221/2同时加快1221/2这个结果看起来很简单但绝大多数刚开始接触的人都会猜错瓶颈是最慢的单台机器所以应该先改进最慢的机器。在这个例子里确实恰好是B但更一般的系统里瓶颈是平均权值最大的环它可能横跨多台机器是一个循环结构而不是单台设备。比如某条产线的瓶颈环可能是机器1加工→机器3加工→回机器1这时候单独优化机器2再猛都没有用。4.4 定时事件图Petri网与极大极小代数的汇合原书在这一章还会讲定时事件图timed event graph。它是Petri网的一个子类每个库所恰好有一个前继变迁和一个后继变迁不存在选择或冲突同时每个库所带有一个时延。在这种结构下变迁的使能时刻恰好可以用dater方程描述而且方程一定是极大极小线性的。事件图的价值在于把Petri网的直观图结构和极大极小代数的计算能力连在一起。你可以在图上画出一个系统直接读出它的系统矩阵算出特征值定位临界环。我的体会是遇到实际问题先在纸上画事件图再写矩阵比直接对着dater方程凑矩阵要快得多也不容易漏约束。原书把这一节放在特征值后面是有道理的——先理解环和瓶颈再看事件图一切就自然串起来了。5. 啃第四章时我踩过的坑五个容易翻车的地方5.1 坑一矩阵里的空约束写成了0构造系统矩阵时两个事件之间如果没有约束对应位置应该写ε-∞而不是0。我最初图省事把所有空位填0结果矩阵特征值变成一个更小的数怎么也算不对。原因在于0在极大极小代数里是不产生延迟的连接等于凭空加了一条两个事件必须同步的约束把系统不存在的等待强加了进来。判断方法很简单删掉这个位置的值递推关系还成立吗如果成立这里就应该是ε。5.2 坑二把极大极小特征向量当成概率分布学过马尔可夫链的人很容易惯性思维看到特征向量就想到稳态分布、想到归一化。极大极小代数的特征向量完全不是这个东西它只是一组相对时间偏移。给所有分量同时加上任意一个常数相对偏移不变特征值也不变在极大极小意义下没有任何归一化要求。真正有意义的是分量之间的差而不是单个分量的值。5.3 坑三dater函数的下标从0还是从1开始写递推时a_k max(a_{k-1}2, b_{k-1}) 这种式子k从1开始和从0开始边界条件完全不同。我第一次自己推导时把a_0设成0结果前两个零件的时刻全偏了。后来养成的习惯是先把k1、2、3这三项的数值手算一遍再倒回去确认递推公式的下标和初始值。代码里就是把a[1]0单独拎出来初始化而不是写一个通用循环去处理k1的情况。5.4 坑四把随机时间硬套进确定性模型这一整章是确定性的时间模型所有寿命、时延都是一个确定的数。我一开始天真地想加工时间明明是随机的把均值代入不就行了实测表明均值代入极大极小系统算出来的吞吐量往往被系统性高估——因为max的本质是取最晚随机波动带来的最晚时刻会比均值情形更晚。随机处理时间是下一章的内容它要用随机定时自动机、广义半马尔可夫过程这套新工具不是简单替换参数。先把确定性模型吃透再去碰随机版本否则两边都会糊。5.5 坑五忽略缓冲区容量约束很多人建模多机串联时只写单机不可重叠和前序加工完成这两类约束忘了缓冲区容量限制。以4.1节的例子来看如果把a_k ≥ b_{k-1}这条去掉吞吐量其实还是1/3因为瓶颈在B但是在别的参数配置下缓冲区约束可能成为临界环的一部分删掉它算出来的吞吐量就会偏高而且排产方案实际执行时会撞墙。我后来检查模型的固定动作就是逐条列出所有谁等谁的关系尤其是跨机器的等待再对照原图确认没有漏约束。6. 第四章和后续章节怎么接上学习路线的一点建议6.1 下一站是随机定时自动机读完第四章最自然的问题是真实系统的加工时间哪有确定值下一章就会把寿命参数从确定数换成随机变量定时自动机升级成随机定时自动机。到那里dater方程的技巧不能直接沿用因为你没法对随机变量做简单的max递推需要转用随机过程视角比如广义半马尔可夫过程和马尔可夫链。所以第四章的确定性假设不是缺陷而是为随机化铺路的简化——先学会在确定世界里算准节拍再考虑波动。6.2 与Petri网、监督控制章节的接口定时事件图是通往更一般定时Petri网的桥梁一旦Petri网里出现选择或冲突dater方程就不再是单纯的极大极小线性形式因为冲突引入了路径不唯一的复杂性。另一方面监督控制理论讨论的是如何施加约束让系统满足规范而时间类规范比如两个事件间隔不得超过10在自动机层面的表示正是以定时自动机为基础的。所以第四章不是一座孤岛它同时给随机系统和控制系统供料。6.3 动手验证的两个具体建议第一一定要把4.1节那种小系统亲手算一遍手写前三次递推再用代码跑前20个周期对照特征值是否稳定在λ。第二自己改一个参数再算一遍比如把缓冲区容量改成2或改成0观察矩阵哪一项变化、特征值变不变。我当初就是把这两件事做完才敢说自己真正看懂了第四章。别人讲十遍不如自己错一遍再改对一遍。另外分享一个小经验啃这章时把事件时间表和图两种视角来回切换。dater方程算得快但没有画面感带权图直观但手工算大系统又很慢。两者互为印证是理解极大极小代数最好的组合。等你能在一张复杂事件图上直接指出临界环的位置第四章就算真正过关了。