熊 濤，丁建完，陳立平

（華中科技大學(xué)國(guó)家CAD支撐軟件工程技術(shù)研究中心，湖北 武漢 430074）

面向半物理仿真的陳述式模型求解方法研究*

熊 濤，丁建完，陳立平

（華中科技大學(xué)國(guó)家CAD支撐軟件工程技術(shù)研究中心，湖北 武漢 430074）

通過(guò)符號(hào)操作和數(shù)值計(jì)算相結(jié)合，提出了一種求解半物理仿真模型的新方法。為了滿足半物理仿真對(duì)實(shí)時(shí)性的要求，在模型編譯階段將代表數(shù)值積分的隱式離散公式插入到仿真模型中，增廣后的方程系統(tǒng)伴隨著非線性方程的出現(xiàn)，需要在積分的每一步對(duì)這些非線性方程進(jìn)行迭代求解，而求解非線性方程的時(shí)間復(fù)雜度隨維度的變大成指數(shù)增加，因此引入代數(shù)環(huán)撕裂減小代數(shù)方程塊耦合變量數(shù)，以滿足實(shí)時(shí)求解對(duì)粒度的要求。最后通過(guò)實(shí)例對(duì)文中提出的方法進(jìn)行了驗(yàn)證。

半物理仿真；微分代數(shù)離散；代數(shù)環(huán)撕裂；實(shí)時(shí)積分

1 引言

仿真技術(shù)是以相似原理、信息技術(shù)、系統(tǒng)技術(shù)及其應(yīng)用領(lǐng)域有關(guān)的專(zhuān)業(yè)技術(shù)為基礎(chǔ)，以計(jì)算機(jī)和各種物理設(shè)備為工具，利用系統(tǒng)模型對(duì)實(shí)際的或設(shè)想的系統(tǒng)進(jìn)行試驗(yàn)研究的一門(mén)綜合性技術(shù)。仿真技術(shù)應(yīng)用廣泛，從早期的航空、航天、火力發(fā)電和核動(dòng)力發(fā)電部門(mén)擴(kuò)展到今天的軍事、電子、通信、交通、冶金、建筑、氣象、機(jī)械、輕工等多種行業(yè)和部門(mén)，其應(yīng)用已滲透到系統(tǒng)生命周期的全過(guò)程［1］。半物理仿真作為工程領(lǐng)域內(nèi)一種應(yīng)用廣泛的仿真技術(shù)，涉及機(jī)電液控多個(gè)領(lǐng)域和學(xué)科?！鞍胛锢矸抡妗庇址Q硬件在環(huán)仿真或物理—數(shù)學(xué)仿真［2］，是指針對(duì)仿真研究?jī)?nèi)容，將控制系統(tǒng)、被控制對(duì)象或系統(tǒng)部件用數(shù)學(xué)模型來(lái)描述，并轉(zhuǎn)化為高速運(yùn)行的實(shí)時(shí)仿真模型，然后與實(shí)物連接成為一個(gè)集計(jì)算機(jī)技術(shù)、數(shù)學(xué)建模技術(shù)及數(shù)據(jù)采集技術(shù)于一體的半物理仿真系統(tǒng)。然后借助物理效應(yīng)模型，進(jìn)行數(shù)學(xué)與物理的實(shí)時(shí)聯(lián)合仿真，通過(guò)仿真試驗(yàn)對(duì)控制系統(tǒng)的控制策略、控制功能以及系統(tǒng)可靠性等進(jìn)行測(cè)試和評(píng)估［3］。目前，半物理仿真已被廣泛地應(yīng)用于軍事、汽車(chē)制造、工程機(jī)械和高校實(shí)驗(yàn)室建設(shè)領(lǐng)域。發(fā)展勢(shì)頭迅猛的多領(lǐng)域物理系統(tǒng)實(shí)時(shí)仿真主要應(yīng)用方式就為半物理仿真，如減速箱［4］、柴油機(jī)控制［5］、汽車(chē)［6］、無(wú) 人汽 車(chē)［7］、鋰 電 池［8］和 航 電模擬［9］等 。

本文對(duì)基于Modelica語(yǔ)言面向半物理仿真的陳述式模型實(shí)時(shí)求解算法進(jìn)行研究，組織結(jié)構(gòu)如下，首先對(duì)半物理仿真及其在工程實(shí)際中的應(yīng)用進(jìn)行簡(jiǎn)單介紹；其次對(duì)半物理模型的數(shù)學(xué)表達(dá)—微分代數(shù)方程DAEs（Differential Algebraic Equations）及其常用的求解方法進(jìn)行分析，給出半物理仿真的難點(diǎn)；然后針對(duì)求解中的困難提出在模型編譯階段對(duì)DAEs進(jìn)行離散，并對(duì)離散后能夠縮減問(wèn)題規(guī)模的代數(shù)環(huán)撕裂進(jìn)行討論；最后給出該方法的優(yōu)勢(shì)和展望。

2 微分方程求解方法

實(shí)時(shí)仿真在仿真軟件應(yīng)用領(lǐng)域發(fā)展的目標(biāo)之一是能夠以快速的采樣率實(shí)時(shí)模擬越來(lái)越復(fù)雜的模型，滿足半物理仿真的特性，也就是在很短的時(shí)間間隔內(nèi)對(duì)物理模型的數(shù)學(xué)描述進(jìn)行求解。工程中通常采用DAEs描述物理系統(tǒng)的動(dòng)態(tài)特性和行為，DAEs通常由相互耦合的微分方程和代數(shù)方程混合而成，由于代數(shù)約束的存在，使DAEs的求解與常微分方程O(píng)DEs（Ordinary Differential Equations）產(chǎn)生了極大的不同，隨著一系列特殊求解方法的提出，DAEs在工程實(shí)踐中得到了廣泛應(yīng)用。其中，指標(biāo)［10］是DAEs的一個(gè)重要屬性，文獻(xiàn)［11］對(duì)低指標(biāo)問(wèn)題進(jìn)行了詳細(xì)研究，高指標(biāo)問(wèn)題可以通過(guò)指標(biāo)縮減轉(zhuǎn)化為ODEs，指標(biāo)縮減代表性的方法有：Gear方法［12，13］、Pantelides方法［14］和啞導(dǎo)方法［15，16］等。

微分方 程 的 求 解 方 法 分 為 變 步 長(zhǎng)［17］和 定步長(zhǎng)［18］兩種，其中變步長(zhǎng)方法需要在每個(gè)積分步對(duì)收斂性進(jìn)行判斷，且必須考慮求解精度對(duì)于步長(zhǎng)的影響，效率較低。“硬件在環(huán)”是半物理仿真的顯著特點(diǎn)，由于回路中接入實(shí)物，因此仿真是實(shí)時(shí)進(jìn)行的，即仿真模型的時(shí)間比例尺和自然時(shí)間比例尺相同，同時(shí)受交互硬件固有頻率的影響，使變步長(zhǎng)算法在半物理仿真中的應(yīng)用很困難。實(shí)際操作中半物理仿真往往使用定步長(zhǎng)的方法，定步長(zhǎng)算法又分為顯式和隱式兩種，對(duì)于常用的顯式方法，當(dāng)模型由機(jī)械、電氣、液壓或者熱力組件耦合而成時(shí)，往往兼有快速和慢速子系統(tǒng)，導(dǎo)致仿真變成剛性問(wèn)題，也就是說(shuō)，模型中的時(shí)間常數(shù)具有較大的跨度，較快的時(shí)間常數(shù)決定仿真的計(jì)算量（步長(zhǎng)），這將帶來(lái)過(guò)多的計(jì)算量，降低了求解效率。使用隱式的求解算法雖然解決了數(shù)值穩(wěn)定性問(wèn)題，并且允許使用較大的步長(zhǎng)，但求解精度制約著步長(zhǎng)，而且使用隱式方法意味著由微分方程離散而來(lái)的非線性方程需要在每一積分步迭代求解，實(shí)時(shí)求解規(guī)模為n的非線性方程系統(tǒng)存在一定問(wèn)題：操作的時(shí)間復(fù)雜度是O（n3），并且問(wèn)題達(dá)到收斂所需的迭代次數(shù)可能隨不同的步長(zhǎng)而變動(dòng)。綜上，如何保證半物理仿真的實(shí)時(shí)性便成為一個(gè)重要的課題，除去仿真計(jì)算機(jī)自身的速度和精度，仿真模型的求解算法是影響仿真速度的一大重要因素，本文給出一種符號(hào)操作與數(shù)值計(jì)算相結(jié)合的半物理模型求解方法。該方法的基本思想為，使用代表求解方法的隱式公式在求解之前對(duì)DAEs進(jìn)行離散，并借助相關(guān)符號(hào)操作實(shí)現(xiàn)求解的規(guī)劃與分解，達(dá)到問(wèn)題分而治之的目的。

3 DAEs的離散

在模型編譯階段將代表數(shù)值積分運(yùn)算的隱式歐拉離散公式插入到方程系統(tǒng)中，實(shí)現(xiàn)指標(biāo)約減后DAEs的增廣，從而將連續(xù)問(wèn)題預(yù)先轉(zhuǎn)化為離散問(wèn)題。完整約束的半仿真模型在指標(biāo)約減后具有下面的表現(xiàn)形式：

其中，E為微分方程，Φ為代數(shù)約束方程，x為狀態(tài)變量，y為代數(shù)變量，假設(shè)在一致性初始條件x（t0）= x0，˙x（t0）=x0，y（t0）=y0下，具有唯一且光滑的解［17］。

DAEs離散步驟：

步驟1 根據(jù)是否具有微分形式將變量劃分為集合x(chóng)和集合y，滿足dim（W）=dim（x）+dim （y），W=［x，y］。

步驟2 為保證方程系統(tǒng)的恰定性，取x為已知量，˙x和y為未知量，對(duì)DAEs進(jìn)行強(qiáng)連通分量SCC（Strongly Connected Components）凝聚［19］，以及拓?fù)渑判?，從而可以得到順序求解的方程子集序列?/p>

步驟3 就xi而言，若˙xi可通過(guò)匹配的方程顯式地求解出來(lái)，則向DAEs中添加如下的離散方程：

其中res（xi）表示將xi當(dāng)作撕裂變量，以迭代的方式求解xi。

若˙xi無(wú)法顯式求解，向DAEs中添加：

步驟4 若與˙xj或者yk匹配的賦值形式方程en處于代數(shù)環(huán)中，即結(jié)構(gòu)關(guān)聯(lián)矩陣中所處的對(duì)角塊維度超過(guò)1，則在en的表達(dá)式中添加一個(gè)新項(xiàng)res（˙xi）或者res（yk），標(biāo)記為撕裂變量，即：

步驟5 如果res（xi）、res（˙xj）或是res（yk）所處的代數(shù)環(huán)不能被撕裂，則移除相應(yīng)的res（）操作符。

步驟6 以x、˙x和y為未知量，重新對(duì)離散后的DAEs進(jìn)行SCC凝聚和拓?fù)渑判颍瑥亩梢缘玫巾樞蚯蠼獾姆匠套蛹蛄小?/p>

下面通過(guò)一個(gè)簡(jiǎn)單的例子進(jìn)行說(shuō)明。

其中，ω、D和k為常數(shù)，f和g為函數(shù)，x、x1和x2為 微 分 變 量，u和y為 代 數(shù) 變 量。若 使 用DASSL［17］直接求解，則需對(duì)全部的5個(gè)變量進(jìn)行迭代。DAEs離散方法處理如下：由于微分變量˙x、˙x1和˙x2均具有顯式表達(dá)，所以使用下面的離散公式（8）～公式（10）對(duì)式（7）進(jìn)行增廣：

增廣后的DAEs如下：

就式（11）而言，為保證方程的恰定性需同時(shí)將˙x、˙x1、˙x2、x、x1和x2視為未知量，消去增廣式中的微分變量˙x、˙x1和˙x2，得：

式（12）的二部圖最大匹配及SCC凝聚如圖1和圖2所示，圖1中雙向邊代表最大匹配，其中需要注意的是式f4中的x1，可以通過(guò)式f3轉(zhuǎn)換為關(guān)于x2的表達(dá)。

Figure 1 Direct constraint graphˉG of equ.（12）圖1 式（12）的約束有向圖ˉG

Figure 2 Cohesion graphˉC of SCCˉG圖2 ˉG強(qiáng)連通分量凝聚圖ˉC

圖2中，強(qiáng)連通分量C1=｛f1|x｝、C2=｛f2|y｝、C3=｛f3|x1｝、C4=｛f4|x2｝和C5=｛f5| u｝，這5個(gè)強(qiáng)連通分量構(gòu)成了代數(shù)環(huán)，當(dāng)選擇變量x作為撕裂變量，反映在圖中為切斷強(qiáng)連通分量C1和C2之間的有向弧，可得順序求解的方程序列C2→C4→C3→C5→C1，方程求解順序與圖ˉC中有向弧的方向相反。

觀察式（13）可知，只需要對(duì)變量x進(jìn)行迭代，即可實(shí)現(xiàn)方程的求解，所以移除標(biāo)記res（x1）和res（x2）。同理，可選擇變量x1或x2作為撕裂變量。下面介紹對(duì)于縮減問(wèn)題規(guī)模具有重要意義的代數(shù)環(huán)撕裂。

4 代數(shù)環(huán)撕裂

為了減小DAEs離散后非線性求解塊的粒度，使用代數(shù)環(huán)撕裂作為后處理，以提高求解效率。代數(shù)環(huán)撕裂的基本思想：利用代數(shù)環(huán)結(jié)構(gòu)上的稀疏性，通過(guò)切斷某些耦合關(guān)系，將某些方程轉(zhuǎn)化為賦值等式的形式，從而縮減需要聯(lián)立求解部分的規(guī)模，提高求解效率?？紤]一般形式的方程組，如式（14）所示：

式（14）為三個(gè)相互耦合的方程，約束表示圖如圖3a所示，這三個(gè)方程組成了一個(gè)強(qiáng)連通分量。撕裂該代數(shù)環(huán)的符號(hào)操作為：打斷組成強(qiáng)連通分量中的某條邊，檢查圖中是否存在新的強(qiáng)連通分量，若存在，繼續(xù)上述的操作；否則，結(jié)束。需要注意的是，每打斷約束表示圖中的一條邊，為保證方程系統(tǒng)的恰定性，需要引入一個(gè)輔助方程和輔助變量，如圖3b所示。

Figure 3 Direct graph of equ.（14）圖3 式（14）的有向表示圖

依次打斷邊（f1，x1）和（f2，x2），得到圖3b所示的形式，方程f'1和變量x'1是打斷邊（f1，x1）時(shí)引入的輔助方程和輔助變量，f'2和x'2是打斷邊（f2，x2）時(shí)引入的輔助方程和輔助變量，符號(hào)操作后的式（14）可以通過(guò)下面的計(jì)算過(guò)程求得：

BEGIN

LOOP

UNTIL|x'1—x1|＜εAND|x'2—x2|＜ε

END

fi_calc_xi表示通過(guò)方程fi計(jì)算變量xi，ε表示迭代收斂所允許的誤差。

下面給出代數(shù)環(huán)撕裂的一般方法，假設(shè)代數(shù)環(huán)具有式（15）指定的形式：

撕裂就是從代數(shù)環(huán)中分離出一部分方程和變量，即｛g1|g1∈g｝和｛z1|z1∈z｝。假定z1為已知量，可由剩余的方程｛g2|g2∪g1=g｝，依次求解剩余的變量｛z2|z2∪z1=z｝?？梢酝ㄟ^(guò)牛頓法求解分裂后的方程系統(tǒng)，該過(guò)程描述如下：

步驟1 選擇合適的撕裂變量z1，并賦初值z(mì)1（0）；

步驟2 通過(guò)剩余的方程g2求出z2，即z2= g2（z1）；

步驟3 計(jì)算殘差δ=g1（z1，z2）；

步驟4 若δ小于收斂誤差ε，結(jié)束，否則，轉(zhuǎn)步驟2。

通過(guò)代數(shù)環(huán)撕裂，可以將耦合方程的維數(shù)從dim（h1）+dim（h2）減小到dim（h1）。然而就撕裂過(guò)程而言，選出合適的z1是十分困難的，不同的選取方式會(huì)帶來(lái)計(jì)算量和誤差上的差異，同時(shí)在實(shí)際操作中必須考慮初值z(mì)1（0）以及h1（z1，z2）迭代的收斂性。為此，采用的方法是取DAEs離散后標(biāo)記為res的變量作為撕裂變量，來(lái)驗(yàn)證能否達(dá)到代數(shù)環(huán)撕裂的目的。

5 實(shí)例分析

本文提出的方法在多領(lǐng)域建模與仿真平臺(tái)MWorks［20］中進(jìn)行實(shí)驗(yàn)，下面通過(guò)笛卡爾坐標(biāo)系下三自由度挖掘機(jī)的轉(zhuǎn)臂模型進(jìn)行驗(yàn)證。

轉(zhuǎn)臂由包括R1、R2、R3、R6、R7、R8、R9、R11、R12、R13、R14、R15和P4、P5、P10在內(nèi)的12個(gè)轉(zhuǎn)動(dòng)副和3個(gè)移動(dòng)副組成，為使存在環(huán)路的模型可解，取R6、R8、R11、R14為切割鉸。記r表示R. phi，用以描述轉(zhuǎn)動(dòng)副的相位轉(zhuǎn)角；p表示P.s，描述移動(dòng)副的伸張長(zhǎng)度。指標(biāo)約減后選擇r7、r12、r15作為狀態(tài)變量的DAEs具有下面的表現(xiàn)形式：

Figure 4 Excavator jib model with three degrees of freedom圖4 三自由度挖掘機(jī)的轉(zhuǎn)臂模型

Figure 5 Directed constraint graph of augmented DAEs圖5 轉(zhuǎn)臂模型增廣DAEs的有向約束圖h

Figure 6 Cohesion graphˉChof SCCh圖6的強(qiáng)連通分量凝聚圖

由圖6可得各強(qiáng)連通分量中需要聯(lián)立求解的變量數(shù)分別為5、3和3，求解順序?yàn)椤?，將各?qiáng)連通分量中標(biāo)記為res的變量作為撕裂變量執(zhí)行代數(shù)環(huán)撕裂，均可以獲得迭代變量數(shù)為1的求解序列，如下所示：

使用隱式歐拉法對(duì)該模型進(jìn)行求解，對(duì)比一般方法（需對(duì)11個(gè)變量同時(shí)進(jìn)行迭代）和文中提出的方法，仿真時(shí)間均為10 s，表1給出了不同方式下積分時(shí)間的對(duì)比。

Table 1 Comparison of different methods in solving time表1 不同方式求解時(shí)間對(duì)比

采用DAEs離散和代數(shù)環(huán)撕裂的方法能夠滿足實(shí)時(shí)性的要求（7.35 s＜10 s），同時(shí)為了保證提出方法的正確性，圖7給出r2的仿真結(jié)果對(duì)比。

Figure 7 Simulation results of r2in two ways圖7 兩種方式下r2仿真結(jié)果對(duì)比

就轉(zhuǎn)臂模型而言，借助DAEs離散和代數(shù)環(huán)撕裂將需要聯(lián)立求解的模型（迭代變量數(shù)為11）轉(zhuǎn)化為三個(gè)代數(shù)方程塊的順序求解，且每個(gè)方程塊的迭代變量數(shù)為1，并且數(shù)值積分時(shí)間小于模型仿真時(shí)間，能夠滿足半物理仿真的要求。

6 結(jié)束語(yǔ)

本文主張?jiān)谀Ｐ途幾g階段采用隱式歐拉法對(duì)DAEs進(jìn)行離散，并通過(guò)代數(shù)環(huán)撕裂來(lái)減小增廣后的方程系統(tǒng)規(guī)模。絕大多數(shù)情況下本文提出的方法能夠?qū)⒄麄€(gè)方程系統(tǒng)分解為多個(gè)可單獨(dú)迭代求解的子系統(tǒng)，從而提高求解效率，滿足半物理仿真的要求，并通過(guò)實(shí)例分析驗(yàn)證了該方法的正確性?？偟膩?lái)說(shuō)，只有在離散后的代數(shù)環(huán)無(wú)法撕裂的情況下，所需迭代的變量數(shù)才與通用積分方法相等。在未來(lái)的研究中希望能夠?qū)Ψ抡婺Ｐ椭械目焖贍顟B(tài)與慢速狀態(tài)進(jìn)行區(qū)分，從而能夠在DAEs的離散中使用顯式的離散公式。

［1］ Liu Rui-ye.Elements of computer simulation［M］.Beijing：Publishing House of Electronics Industry，2011.（in Chinese）

［2］ Aranya C，Xin Yu-feng.Hardware-in-the-loop simulations and verifications of smart power systems over an exo-GENI testbed［C］∥Proc of the 2nd GENI Research and Educational Experiment Workshop，2013：16-19.

［3］ Kang Feng-ju，Yang Hui-zhen，Gao Li-e.Technology and application of modern simulation［M］.Beijing：National Defense Industry Press，2006.（in Chinese）

［4］ Casella F，Donida F，Lovera M.Beyond simulation：Computer aided control system design using equation-based object oriented modeling for the next decade［C］∥Proc of the 2nd International Workshop on Equation-Based Object-Oriented Languages and Tools，2008：35-45.

［5］ Song Bai-ling.Semi-physical simulation technology research of diesel electronic control system［D］.Harbin：Harbin Engineering University，2009.（in Chinese）

［6］ Mattsson S E，Elmqvist H，Otter M.Physical system modeling with Modelica［J］.Control Engineering Practice，1998，6 （4）：501-510.

［7］ Tang Guo-ming.Design of hardware-in-the-loop simulation system for driverless vehicle［D］.Hefei：University of Science and Technology of China，2012.（in Chinese）

［8］ Luciano S，Ines C，Manuela G.A design methodology for semi-physical fuzzy models applied to the dynamic characterization of LiFePO4 batteries［J］.Applied Soft Computing，2014，14（B）：269-288.

［9］ Jiang Zhan-si.Parameter optimization and inference solving for models based on Modelica［D］.Wuhan：Huazhong University of Science&Technology，2008.（in Chinese）

［10］ Unger J，Kr?ner A，Marquardt W.Structural analysis of differential-algebraic equation systems-theory and applications ［J］.Computers&Chemical Engineering，1995，19（8）：867-882.

［11］ Negrut D，Jay L O，Khude N.A discussion of low order numerical integration formulas for rigid and flexible multibody dynamics［J］.ASME Journal of Computational and Nonlinear Dynamics，2007，4（2）：1-11.

［12］ Gear C W，Petzold L R.ODE methods for the solution of differential/algebraic systems［J］.SIAM Journal on Numerical Analysis，1982，21（4）：716-728.

［13］ Gear C W.Differential-algebraic equation index transformations［J］.SIAM Journal on Scientific and Statistical Computing，1988，9（1）：39-47.

［14］ Pantelides C C.The consistent initialization of differentialalgebraic systems［J］.SIAM Journal on Scientific and Statistical Computing，1988，9（2）：213-231.

［15］ Mattsson S E，S?derlind G.Index reduction in differentialalgebraic equations using dummy derivatives［J］.SIAM Journal on Scientific Computing，1993，14（3）：677-692.

［16］ Feehery W F，Barton P I.A differentiation-based approach to dynamic simulation and optimization with high-index differential-algebraic equations［C］∥Proc of the 2nd International Workshop on Computational Differentiation，1996：239-252.

［17］ Brenan K E，Campbell S L，Petzold L R.Numerical solution of initial-value problems in differential-algebraic equations ［M］.Philadelphia：SIAM，1996.

［18］ Mathews J H，F(xiàn)ink K K.Numerical methods using Matlab ［M］.London：Prentice Hall，2004.

［19］ Xavier A.On the complexity of strongly connected components in directed hypergraphs［J］.Algorithmica，2014，69 （2）：335-369.

［20］ http：∥www.tongyuan.cc/.

附中文參考文獻(xiàn)：

［1］ 劉瑞葉.計(jì)算機(jī)仿真技術(shù)基礎(chǔ)［M］.北京：電子工業(yè)出版社，2011.

［3］ 康鳳舉，楊惠珍，高立娥.現(xiàn)代仿真技術(shù)與應(yīng)用［M］.北京：國(guó)防工業(yè)出版社，2006.

［5］ 宋百齡.柴油機(jī)控制系統(tǒng)半物理仿真技術(shù)研究［D］.哈爾濱：哈爾濱工業(yè)大學(xué)，2009.

［7］ 唐國(guó)明.無(wú)人駕駛汽車(chē)半物理仿真系統(tǒng)設(shè)計(jì)［D］.合肥：中國(guó)科學(xué)技術(shù)大學(xué)，2012.

［9］ 蔣占四.Modelica仿真模型的參數(shù)優(yōu)化及推理求解研究［D］.武漢：華中科技大學(xué)，2008.

熊濤（1990），男，河南確山人，博士生，研究方向?yàn)槎囿w系統(tǒng)動(dòng)力學(xué)、多領(lǐng)域建模與求解。E-mail：xiongtao39@163.com

XIONG Tao，born in 1990，Ph D candidate，his research interests include dynamics of multibody systems，multi-domain modeling and solving.

A method for solving semi-physical simulation declarative models

XIONG Tao，DING Jian-wan，CHEN Li-ping
（National CAD Support Software Engineering Research Center，Huazhong University of Science and Technology，Wuhan 430074，China）

We propose a new method for solving semi-physical simulation models using a mixed symbolic and numeric approach.In order to meet the real-time requirement of semi-physical simulations，implicit discretization formulae representing the numerical integration algorithm are inserted into the DAEs symbolically at compile stage.Then nonlinear equations will appear in the augmented equation system and all these nonlinear equations should be solved together at each integration step.In order to meet the fine granularity required for solving the models in real-time，tearing algebraic loop is introduced.After that the dimensions of nonlinear equation blocks can be reduced as the time complexity of calculating nonlinear equations increases exponentially with the growth of dimensions.Finally，an example is given to show that the proposed method is not only easy to implement but also efficient.

semi-physical simulation；discrete DAEs；tearing algebraic loops；real-time integration

TP301

A

10.3969/j.issn.1007-130X.2015.08.018

1007-130X（2015）08-1540-06

2014-06-27；

2014-10-11

國(guó)家科技支撐計(jì)劃資助項(xiàng)目（2012BAF16G02）

通信地址：430074湖北省武漢市華中科技大學(xué)國(guó)家CAD支撐軟件工程技術(shù)研究中心

Address：National CAD Support Software Engineering Research Center，Huazhong University of Science and Technology，Wuhan

430074，Hubei，P.R.China