仿真系统中DAE 求解技术现状

杨文强，吴文渊，陈经纬，冯勇

（1.中国科学院重庆绿色智能技术研究院，重庆 400714；2.中国科学院大学 重庆学院，重庆 400714）

企业为了提高在全球化市场中的核心竞争力，提高产品的开发效率，降低产品的开发成本，优化产品性能和结构设计，计算机建模仿真技术将逐渐成为其必不可少的工具[1]。基于对仿真技术的需求，应运而生了数字化设计与制造技术，包含了基于云计算和互联网的用户需求挖掘与新产品开发策划、计算机辅助设计（Computer Aided Design，CAD）、计算机辅助工艺设计（Computer Aided Process Planning，CAPP）、计算机辅助工程（Computer Aidedengineering，CAE）、计算机辅助制造（Computer Aided Manufacture，CAM）、数字样机（Digital Mock-Up，DMU）等[2]，并结合产品“V”模式开发，能够大幅缩短产品研发和测试的时间，有效避免不合理的设计造成的成本浪费。

近年来，随着计算机仿真技术的发展，通过建立数学模型和软件算法相结合的形式来解决工程问题的方法得到了广泛应用。实现从设计到生产的过渡过程中软、硬件结合的桥梁——硬件在环技术（Hardware in the Loop Simulation，HILs）逐渐变得至关重要。它利用计算机建模仿真技术对受控对象的运行状态进行数字孪生，并依靠处理器实时输出仿真结果，通过I/O 接口与控制器连接，对控制器进行统的测试，有效地提高了系统的安全性、可行性，降低了长期设计成本[3]。新技术的提出，使仿真技术逐步发展出了离线（Off-line）仿真、在线（On-line）仿真及离线和在线相结合的仿真模式，极大地促进和保证了虚拟还原现实的可靠性和准确性。

1 仿真建模中的DAE 方程

当前，伴随着大数据的发展，作为世界制造工厂的我国，正面临着数字化、智能化改革与创新的挑战和机遇。新产品日新月异，更加广泛地面向个性化定制的需求，在产品结构和功能上更加复杂和丰富，诸如机器人、无人汽车等多智能体形式的现代高科技产品，由于其鲁棒性、可靠性、高效性、可扩展性等特性，在计算机网络、机器人、制造业、交通控制、虛拟现实等方面得到广泛应用。与传统的机电系统、电气系统等不同，现代高科技产品是更高层面上的多个学科领域子系统集成于一体的复杂大系统，即多领域（multi-domain）耦合的复杂大系统。该系统具有多体系统、多物理场、多学科交叉融合的特征，是现阶段制造技术发展的必由之路。要实现对现代高科技产品的设计，必然离不开利用计算机建模仿真技术对其进行多领域的协同仿真，以实现产品整体性能的高精准表述[1]。

1.1 仿真建模软件的分类

目前，市场上仿真软件主要分为专业的仿真软件和通用的仿真软件。专业的仿真软件，如多体机械系统仿真的ADAMS、电磁场仿真的ANSOFT 和流体力学仿真的ANSYS 等[1]，由于其专注性太强，缺乏其他领域的可扩展性，因此不适合应用在多领域统一建模中。具有代表性的多领域统一建模通用仿真软件主要有 MatlabSimulink，以及基于 Modelica 语言的SimulationX、MapleMaplesim、System Modeler、Dymola、Openmodelica、Mworks 等。其中，以MatlabSimulink软件为代表的因果建模需要明确定义元件的输入和输出，而以开源的Modelica 语言软件为代表的非因果建模则无需考虑元件之间的信号流。同时，MatlabSimulink 因其在矩阵计算中的突出特点，使用最为广泛，同时兼具了计算能力强，计算效率高等特点，但建模过程复杂。基于Modelica 语言的非因果建模方便将模型投影到数学表达式，更加方便使用者直观地理解，并且在符号计算层面上保障了计算的可靠性。特别地，还需要考虑到MatlabSimulink 并非开源软件，Mathworks 公司于2020 年初受美国政府的干预对中国部分高校和科研机构进行了制裁[4]。

1.2 仿真建模求解流程

仿真建模求解流程基于Modelica 语言仿真器的编译求解过程，大致分为编译、分析优化和仿真求解3 个阶段[5]，见图1。

图1 基于Modelica 语言的仿真器求解过程Fig.1 Solving process based on Modelica language simulation

其中，编译阶段可进一步分为词法分析、语法分析、语义分析和平坦化4 个部分，主要完成将Modelica模型的程序语言源代码进行计算机数学表达方程化处理的过程，利用词法规则、语法规则进行程序验证和分析，提取出其中的方程信息特征，平坦化后得到一组方程及其对应的常量、参数和变量。

分析优化阶段是利用数学方法完成对方程求解前的相容性分析和模型简化的预处理工作。其中，相容性分析是利用二部图等方法检查是否完美匹配，即方程和变量关系是否相等，也就是系统是否为恰约束系统。同时，对过约束和欠约束问题，能够协助用户高效地实现问题分析和解决方案拟定，帮助实现模型的校正。模型简化有利于方程的快速求解，涉及方程表达式的符号约简、DAE 指标约简、方程系统的分块化处理，以及转化为状态空间表达式等技术，得到一个DAE 系统的可顺序求解的方程子集序列。

仿真求解阶段与ODE 求解不同，DAE 求解中不可或缺的重要环节是对代数约束方程进行一致性初值（Consistent Initial Value）的求解。求解器根据DAE方程子集求解序列的类型筛选出数值求解包中最优的求解函数，并按照求解函数格式将方程子集求解序列自动生成基于C 语言的求解算法的可执行代码，通过编译器进行执行并输出结果。

Simulink 仿真求解过程与基于Modelica 语言的仿真器的求解过程类似，也可分为仿真建模、编译、分析优化和求解这4 个过程。Simulink 仿真的区别在于：仿真建模阶段，Simulink 不仅需要确定传送宽度、数据类型、计算块参数以及分配内存，还需要明确信号流方向和采样时间；分析优化过程中，为了更好地发挥矩阵计算的效率，需要用户根据情况选择合适的求解器对其进行数值求解；在求解的过程中，计算是离散分块的数据流迭代求解的过程。

目前，编译阶段及分析优化阶段的相容性分析和模型简化都属于计算机程序语言的设计工作，都已经有了很好的解决方案。在指标约简及求解阶段，虽然也有比较成熟的解决方案，但仍然存在理论上的漏洞和巨大的可发展和优化的空间。

目前，在软件的底层算法上，都需要将多领域物理模型平坦化为数学模型的微分代数方程，进而对DAE 进行求解。因此，在数学上，多领域统一建模的仿真问题就是对应DAE 方程的求解问题[5]。

DAE 方程的一般形式可以描述为：

式中：t为自变量，t∈I；I为非空区间，I⊂R；x为n维因变量，x=x(t) = [x1(t), …,x n(t)]；x(k)表示x(t)的k阶导数，1≤k≤l，k∈Z。

与ODE 不同的是，DAE 由含有x的最高次导数x(l)的微分方程和不含x的最高次导数x(l)的代数方程2 个部分构成，故而DAE 的解取决于因变量x和它的导数，而不像ODE 的解那样仅取决于x自身。

DAE 方程的求解实质上是将DAE 方程中关于因变量x导数的隐含约束方程找出来，即进行指标约简使其转化为ODE，然后再对ODE 进行数值求解。其中，涉及指标约简、一致性初始值、常微分方程数值精确求解等几个关键技术。

2 DAE 指标约简

DAE 系统的指标按照不同的定义和用途等分为微分指标、结构指标、Kronecker 指标、Perturbation指标等[6]。虽然它们都在一定程度上代表了DAE 的求解困难和求解精度，指标越高求解越困难，求解精度越差，但究其本质最终目的都是求微分指标。微分指标是用来度量它与ODE 之间的“距离”，也就是将其转换为ODE 系统而必须对其中部分或全部方程求微分的最小次数[7]。理论上，经过恰当次数的微分操作后，DAE 中所有次数导数都可以用新变量替换，得到等价表达的多项式系统[8]。高指标（指标≥2）的DAE 系统求解不可避免会用到数值微分运算，不同于数值积分的情况，步长作为数值微分的分母项，减小步长会导致数值微分误差的急剧增大，经常会引起数值求解的不稳定，造成结果错误，因此无法通过减小步长来实现误差控制[9]。现常用的DAE 求解器通常只能直接求解低指标的问题或特殊形式的高指标问题。

2.1 常见DAE 指标约简方法

对高指标DAE 的数值求解问题，一般有2 种思路：一种是直接对DAE 进行求解，但容易出现降阶、漂移和不稳定的现象；另一种是先将高指标DAE 约简为指标为1 的DAE 或ODE，然后再对约简后的系统方程进行求解[6]。Gear 方法[10]、Pantelides 方法[8]、哑变量方法[11]、负权二部图法[1]、加权二部图法[12]、Pryce[13]方法是目前具有代表性的几种指标约简方法，可对一般形式的DAE 系统实现指标约简。

Gear 方法是一种纯符号的计算方法，反复地对代数方程进行微分及部分低阶微分变量替换，最终将DAE 系统转化为ODE 系统。虽然该方法能完全实现指标约简，但符号算法复杂、计算效率低，且容易出现一些不必要的微分。

Pantelides 方法是一种结构化分析方法，它通过分析方程系统结构对应的雅可比矩阵，寻找出其中奇异的最小结构子集，然后对最小结构子集微分，如此往复直到雅可比矩阵非奇异。该方法由于采用数值方法实现，与Gear 方法的符号方法相比，其计算量大幅降低，且避免了重复性微分，时间复杂度低，目前应用在Dymola 和Maplesim 仿真软件中。

哑变量方法通过对函数变量引进新的哑变量来解决Pantelides 方法中动态变量的选择问题，本质上是Pantelides 方法的一种改进版。

负权二部图法和加权二部图法是Pantelides 方法的补充方法，2 种方法相似，基于赋权二部图的直接实现高指标DAE 的约简，且与Pantelides 方法具有相同的时间复杂度。

Pryce 方法是一种Pantelides 方法的等效方法，基于DAE 系统符号矩阵（Signature Matrix）的指标约简方法，通过求解对偶问题的整数规划问题（Integer Linear Programming，IPL），能够直接地获取系统的最优偏移量和筛选出所需指标约简的方程。

其中：DAE 系统F对应的n×n维的符号矩阵σ中的元素构造如下[13]：

对偶问题的整数规划问题的形式如下：

式中：c为方程偏移量，c=[c1,…,cn]；d为导数偏移量，d=[d1,…,dn]，δ(F)为最优解问题的最优值。

2.2 DAE 指标约简修正方法

DAE 系统指标约简的成功与否取决于指标约简后的雅可比（Jacobian）矩阵是否奇异，如果不奇异，则指标约简成功；其中，雅可比矩阵为DAE 系统指标约简后因变量最高次导数x（l）的系数矩阵。值得注意的是，将DAE 系统方程经过非奇异的线性矩阵A进行线性重组处理，即A·F，或对其中某些方程进行平方处理，即=0，形成新的DAE 系统，采用上述基于结构的指标约简方法，如Pantelides 方法、Pryce 方法等，很可能会出现指标约简不成功的情况，究其本质，也是由雅可比矩阵奇异导致的。对带平方的DAE 进行指标约简后，生成的ODE 系统中会存在冗余方程，影响数值求解的进行。

雅可比矩阵奇异的情形可以分为2 类：一类是变量在求解区间的某些孤立的点上的奇异，但整个区间上是连续可积的，如例1 所示，这一类的奇异不会影响DAE 系统求解的质量，通常在求解过程中可以通过坐标系变换等方法来摆脱这一类奇异的影响；另一类则是退化，即在整个求解区间上都是奇异的。进一步，退化问题又可以按照能否通过符号计算的方式约简为0，分为符号退化（如例2 所示）和数值退化（如例3所示）2 类，其中，符号退化问题通常表现为结构上可约导致的雅可比矩阵奇异，包含结构奇异（Structurally Singularity）和完全奇异（Identically Singularity）。

例1，活塞缸的曲柄滑块机构（见图2）如图2所示，根据欧拉-拉格朗日方程分析得到该曲柄滑块机构的DAE 方程如下。

图2 曲柄滑块机构Fig.2 Slider-crank mechanism

假设转动惯量J1= 1、J2= 2，杆长l1= 1、l2= 2；当燃料燃烧对活塞（滑块）做功，推动活塞运动，其运动规律可控，满足δ˙˙ - sint=0，此时有该DAE 的雅可比矩阵行列式值为 4cosθ1·sinθ2·cosθ2+ 4sinθ1·(cosθ2)2，当2 根杆重叠时，曲柄滑块机构处于死点位置，雅可比矩阵奇异。此位置如果为初始启动位置，势必造成机构无法启动，仿真求解失败，可以采取额外的措施增加方程数量修正雅可比矩阵的奇异；当运动过程中经过该位置，则由惯性会继续运动，避免卡死，对位置的仿真通常会由计算误差或步长而忽略不计。

例2，符号退化问题可以分为以下3 类。

1） 结构奇异（不完美匹配）[14]

2） 结构奇异（完美匹配）[15]

3） 完全奇异[14]——非线性单摆模型

对二部图不具有完美匹配的结构奇异问题，如例2 中1）所述，已经不是一个单独的DAE 系统了，多个系统之间相互独立，应对其单独进行仿真建模，而无共同仿真的实际意义。对二部图具有完美匹配的结构奇异问题，如例2 中2）所述，可以通过符号约简的方法对其进行预处理，消除其中的冗余部分[16]。

对于完全奇异问题，如例2 中3）所述，虽然雅可比矩阵J奇异，但仍然存在指标约简方法未能找出全部的隐含约束条件的情况，因此需要进一步对其进行修正。最早的修正方法是Murota[17]提出的通过组合松弛算法及其改良方法[18-19]，对线性的DAE 系统进行重组表达。随后，在其修正框架上发展出改进的修正方法：线性组合方法（LinearCombination，LC）、表达式替换法（ExpressionSubstitution，ES）[15]、替换法（Substitution Method）和增广法（Augmentation Method）[14]。

基于组合松弛方法中提出的修正框架概括为以下3 步。

1）计算等效整数规划问题的最优解(c,d)，如果该问题没有解，则DAE 将无法实现完美匹配，算法将以失败告终。

2）判断雅可比矩阵是否奇异，如果不是，则返回ODE 系统F(c)。

大多数的退化问题都是符号退化的情形，只有少数，例如在含参数的微分代数方程中，参数在取值范围变化会引起雅可比矩阵的变化，当取到某些特殊参数时，难以避免会出现数值退化问题，此时，现有的修正方法将不能胜任。

图3 共线弯矩作用下梁的叠加变形Fig.3 Superposition deformation of beams under collinear bending moments

2.3 一致性初值问题

将DAE 系统进行指标约简后，只是完成了对隐含约束方程部分的处理，但还缺乏对初值问题或边界值问题相关约束条件的描述，因此还需要解决一致性初值问题，才能将DAE 系统完全地转换成为一个ODE 系统，同时，结合例3，可知初值问题的求解有利于数值退化问题的分析。

在仿真建模中，系统的初始状态通常可以不同，如曲柄滑块机构在正反方向上都可以进行启动，但其运动的状态——拉伸或压缩，却是截然不同的。如果只是对一个初始方向进行分析和优化，则所得的结果都是局部的而非全局的。

例如，对代数方程形如双曲线方程x2-y2-1 =0 这一类，如果采用牛顿迭代等数值方法对其进行求解，势必会造成部分一致性初值信息的丢失。而这些信息有可能会关系雅可比矩阵的奇异性，如果单纯的放任不管，则可能导致DAE 系统与一致性初值的不匹配，造成求解失败。

针对类似的多项式方程求解问题，常采用Sommese、Wampler、Versheclde 等提出的数值代数几何学[20-21]的同伦延拓求解方法（Homotopy Continuation Mothed）来进行处理。对于非线性较强（如含三角函数）的情况，该方法仍然无法实现求解。

3 ODE 求解器

在机械动力学和控制工程等领域的模拟仿真中，系统在某时刻后的状态变化，通常被描述为线性常微分方程的初值问题。如果常微分方程在时间t上连续，在x上Lipschitz 连续，那么它在边界条件的领域范围内一定存在唯一的精确解*x。目前，在线性常微分方程的初值问题和[22]边界值问题[23-24]上，已有非常成熟的数值求解方法和p阶误差的分析方法，并且开发了很多面向应用的求解器，例如Matlab[25]和Maple[26]等。据了解，关于如何严格给出由求解器数值解插值得到的近似解与精确解的距离上界，以及如何在数值解的基础上提高近似解的质量等方面的工作还比较少。

3.1 ODE 误差估计方法

对误差估计而言，通过残差对其进行界定具有许多优点，最大的优点就是它能够用来评估整个求解区间内的误差，而不只是每个网格格点上的局部误差[27]。研究人员通常喜欢将求得的数值解代回原方程求出残差[28-30]，再进行误差估计。Constantinescu[31]在文章中提出利用时间步进策略结合残差给出一个p阶“全局误差”的估计方法，能够更好地提高误差估计的准确性，但本质上只是在格点上的误差估计，而且舍去了高阶项，不能说是准确的误差估计。Enright[32]在文章中提出了合理的估计残差范数方法，但没有考虑每个点导数的数值误差，这些导数的不准确进一步导致残差定义的不准确，从而导致误差估计的不准确。另一些文献[28,30,33]在考虑导数问题的误差估计方面做了大量工作，但其中大部分都在努力解决局部误差问题。

一般情况下，人们是无法求出常微分方程的精确解的，但幸运的是，可以通过很多数值方法求解器求得数值解，如Matlab 等，可以进一步减小步长来逼近精确解，基于此，很多人在这方面做了重要的研究工作[34-36]。无论步长取得多小，局部误差控制得多好，人们仍然无法回答数值解和精确解之间到底相差多少这个问题。全局误差估计方法[37-40]能够部分地回答这个问题，但是这些方法计算量太大且条件过于苛刻。文献[41]中提出利用Lipschitz 常数来协助误差上界的估计，给出了误差估计上界：

式中：h为步长；M为的估计上界；L为Lipschitz 常数，由于Lipschitz 常数为正，导致估计上界随时间增长呈指数增长，且增长速率过大，不具有实用性。文献[29]通过推导条件数和残差的关系来估计全局误差，但求解条件数需要对齐次常微分方程的基本解（Fundamental Solution）进行求解，在多数情况下这一点几乎无法实现。专门针对稳定系统的常微分方程，文献[42]提出利用Lyapunov 定理构造全局误差控制方程来限制误差上限，能够起到误差估计和提高解的精度的作用，但构造全局误差控制方程也需要基本解的求解。文献[43]中，提出了一种区间求解的方法，仅适用于初值问题求解，该方法基于泰勒级数展开，估计的范围和计算精度取决于级数展开的项数，精度提高会造成计算量急剧增大；同时，经实验验证，区间求解方法需要大量的预处理时间；且在相同的步长和阶数的情况下，区间求解方法的精度要低于龙格库塔等数值求解算法。文献[44-45]利用残差和精确解之间的关系，通过理论推导和证明构造了误差估计上界，能够准确地实现线性ODE 初值问题的误差上界的估计，避免了Lipschitz 常数的指数增长问题和求解基础解的问题，但针对非线性和边界值问题还缺少进一步的研究。

在数值求解中，解的精度通常取决于步长，但更小的步长意味着更大的计算量，如何在固定步长的情况下减小全局误差是一个很有难度的问题，目前还没有太多相关的研究。“缺陷修正”（Defect Correction）方法[46-47]，类似于牛顿迭代方法求解非线性方程组，将残差代回方程进行迭代数值求解，可用于非线性常微分方程，但它们都只是考虑了在网格点上的局部误差。同时，在计算效率上，虽然现有的数值求解方法都很快，但仍然可以考虑利用技术手段将计算效率进一步提升。

3.2 ODE 数值求解方法

现有的ODE 系统的求解器有很多，大致可以分为3 类：数值求解器、符号求解器和区间求解器。数值求解器主要基于欧拉法及其改进格式、龙格库塔格式、亚当斯格式等[48]，其优势在于计算效率高，目前应用最为广泛，常见的为Matlab 提供的求解器，针对不同类型的方程与不同问题，比较典型的有针对常微分方程的ode45、针对时滞微分方程的dde23、针对边值问题的bvp4c 和bvp5c 等。符号求解器是基于微分方程通解公式，主要应用在Malpe 的求解器中，比较典型的求解器有dsolve、linearsol 等。目前只能针对部分格式的常微分方程进行求解，如齐次方程、线性方程、Bernoulli 方程、Riccati 方程等格式。

Nedialkov 等[43]的区间求解方法需要采用高阶的泰勒级数展开对微分方程进行预处理，然后才能够得到网格点所在的一个时间区间上的解区间，该网格点上的精确解属于该解区间，目前应用在VNODE-LP区间求解器。由于没有将数值解插值作为近似解，因而区间求解方法只能给出网格点所在区间上解得的误差估计区间，无法给出全局的误差估计。

有限差分方法是一种应用范围更广的数值求解方法，该方法不仅能够用来求解一般的ODE、DAE、偏微分方程（Partial Differential Equation，PDE），且不区分初值问题与边界值问题，还能够用来求解ODE或DAE 相关的优化问题。其基本思想是通过对定义域进行网格划分，将原问题在网格点上进行离散化，然后用代数差分方程组去逼近原来的微分方程，同时需满足对应原问题的在离散点处的边界条件，典型的有限差分方法包括配置法和有限元法[49]。

文献[45]利用残差二范数和精确解插值相结合的符号推导方法，给出了一种线性ODE 和部分非线性ODE 的优化求解方法，同样能够解决初值问题和边界值问题的求解。该方法首先需要在符号层面上完成ODE 方程的预处理后才可赋值进行计算。虽然采用离线和在线的方法能够极大地加速特定类别方程的计算效率，但针对非线性强的ODE 方程求解还需要进一步研究。

随着深度学习技术的发展，利用神经网络架构求解常微分方程[50]、偏微分方程[51]成为了可能，这也是目前比较热的研究领域，是一个值得研究的方向。Chen 首先提出了将残差网络连续变化为微分方程的思想，能够高精度地对常微分方程的初值问题进行求解，以及在此基础上演变而来的[52-53]，并非通用的求解模型，针对不同类型的ODE 都需要采用训练集对神经网络进行训练，效率较低。同理，基于数据驱动的回归分析方法[54-55]来求解偏微分方程也存在这样的问题。

4 DAE 系统中的参数设计

当前的仿真建模软件还处于发展的初期阶段，大多数情况下只是考虑如何完善对一般形式的DAE 系统的求解功能，而对含参数的DAE 系统问题涉及较少。诚然，一般形式的DAE 系统的求解功能是含参数的DAE 系统问题求解的基石。对产品进行优化设计、对仪器设备进行优化控制、对原材料和制备过程进行优化调节，是提高行业竞争力的最直接方法。如汽车领域，在设计时可以利用仿真器优化行星齿轮机构特征参数，使之与动力部件的最高转速和整车设计最高车速匹配，使发动机工况得到了改善，整车燃油经济性和排放都得到大幅提升[56]。究其本质，就是对DAE 系统进行参数设计，即参数优化，最常见的还有对动力系统的PID 参数调节。而对产品参数化设计，如果采用传统的“V”模式开发，势必会造成时间和金钱的浪费，因此采用参数优化方法结合仿真建模的方式会是未来仿真建模领域发展的主要趋势。

对含参数DAE 系统而言，目前应用最多的优化方法都是局部优化法，具体包括：利用极大值原理将含参数的DAE 模型一般化，然后再利用DAE 求解器进行求解的变分法；将优化含参数的DAE 模型为非线性规划模型，进而进行优化求解的联立法；有限差分方法[49]。采用局部优化方法收敛得到的局部最优很可能是一种假象，实际情况是人们可能从来没有找到过“局部最优”。

针对含参数DAE 系统的全局最优化方法方面的研究很少，可以参考在含参数ODE 问题求解的确定性方法和随机性方法。随机性方法（如启发式方法），通过不断搜索的方式，逐个比较极值，力求遍历全局，但计算量大、效率低下，典型的有蚁群算法[57-58]、退火算法[59-60]、粒子算法[61]等。这类方法没有严格的理论支撑，能够收敛到极值，但不能保证获得的最优解就是全局最优解。确定性方法是基于函数的解析性质（如凸性、稠密性、单调性和Lipschitz 等），然后利用这些性质确定含参数DAE 系统的界，进一步确定全局最优解[62]，主要有分支定界法、整数全局寻优、聚类和隧道方法等。分支定界算法是将原问题通过松弛技术转化为简单的凸优化问题[63]，已经实现部分含参数的ODE 系统求解中的应用[62-64]，但其在计算过程中需要求解ODE 系统的基本解来确定凸函数的上界和下界，因此只能针对特定的ODE 优化系统进行求解。文献[65]中提出的分支定界方法框架下新的凸松弛方法能够处理多项式形式的广义几何规划等问题，但其仅限于特定形式的多项式问题。目前，关于分支定界方法的研究主要集中在凸松弛和分支策略上，其目的是通过紧度来保证解的质量，通过裁剪分支来保证解的效率；其中，最重要的环节就是凸松弛技术，包括：线性松弛[66]、凸二次松弛[67-68]、拉格朗日松弛[69]、二阶锥松弛[70]，以及半正定松弛[71]等。

5 未来研究方向

多领域统一建模仿真软件的目标就是通过计算机技术实现虚拟现实模拟及数字孪生，协助开发人员便捷、可靠、有效地进行新产品的开发设计。为了更好地实现这一目标，必须要满足DAE 求解的可靠，而根据之前所述内容可知，现有DAE 求解技术存在诸多局限，对仿真软件求解器进行优化，构建一个更加精确和完备的DAE 求解架构，显得非常有必要。未来在这一领域的研究可以归纳为以下几个方向。

1）在现有的仿真模型、编译器和分析优化器等环节引入模型的校验、验证和确认（Verification，Validation and Accreditation，VV&A）技术，提高DAE仿真模型的正确性，解决其中的过约束和欠约束方程问题，这样能够提高DAE 求解的可行性、真实性和可靠性。值得考虑的是：可以采用硬件在环技术、机器学习技术（如对抗生成学习网络、强化学习网络等）来实现；特别是VV&A 技术在针对大型复杂的DAE仿真系统时，在精度和置信度的提高上，可以做进一步的研究。

2）虽然高指标DAE 系统的指标约简技术已经较完善，但在DAE 系统出现退化（尤其是数值退化）的情形时，以及在DAE 中存在方程的线性重组或平方的情形时，现有的指标约简及其修正方法可能会失效，亟需对现有技术进行修正或补充。

3）一致性初值问题是求解DAE 系统的关键问题之一，对于单分支的情况，采用数值求解约束条件方程是完全能够胜任的，但对于多分支的情况，有必要遍历每个分支，考虑DAE 多个状态的初值，数值求解方法只能得到个别的解。目前可以采用的解决方法有同伦连续方法等，但缺乏其在DAE 一致性初值上的相关理论推导和证明，以及在非多项式约束条件方程会失效。

4）针对ODE 数值求解，现有仿真软件中的数值求解器已经非常成熟，但对非线性程度高的大规模复杂模型上还可以进一步研究，可以考虑结合神经网络技术；全局误差可控计算在高精度过程控制中有非常广的应用前景，但对其的研究还很少，目前的研究局限在近似线性ODE 系统，且符号计算不利于大规模求解，短期来看，可考虑利用泰勒级数展开等技术推广到非线性的复杂模型上。

5）在ODE 解的误差估计上，目前在局部误差和全局误差的估计上存在很多有效可行的方法，但全局误差估计方法的研究还局限在线性的初值问题上，可考虑利用线性化技术或初值转换技术推广应用到非线性和边界值问题的情形，这也是一个不错的研究方向。

6）含参数DAE 系统的优化问题任重而道远，特别是全局参数优化问题。短期来看可以利用指标约简技术，将DAE 方程转化为ODE 方程，再利用现有的优化方法（如分支定界法）实现特殊形式的微分代数方程的参数优化（如低指标、线性、常系数等），这类问题在控制工程中有着广泛的应用前景。长期来看，需要在理论上进行创新，如在Singer[64]研究的基础上，避开基础解的求解，从根源上解决这一类问题。

6 结语

仿真建模的关键就是建立DAE 系统和求解DAE系统，前者是产品设计开发人员根据产品的物理、数学和经验等构建的有机结合，随着模型的改变而改变；后者是仿真建模软件的核心，其中涉及指标约简技术、数值求解技术、误差分析技术及优化求解技术等。现有的指标技术和基于此的修正方法能够满足绝大多数DAE 系统转化为ODE 系统的需求，但针对特殊形式（如数值退化问题和多分支一致性初值问题）还需要进一步研究。同时，现有的数值求解技术也相当成熟，用户可以根据需求去合理地选取，但在连续区间误差函数可积的情况下，全局误差最小化也是一个可以研究的方向。其次，现有的局部误差估计方法和全局误差估计方法都比较准确，如何进一步控制误差估计上限，以及针对边界值、非线性问题进行准确误差估计，同样是一个不错的研究方向。最后，含参数的DAE 问题或含参数的ODE 问题的求解仍然存在诸多的困难和挑战，关系到产品个性化发展的未来。