基于Matlab與三維空間幾何變換的車輛動態包絡線獲取方法研究

戶攀攀 劉建新 蔡久鳳 劉志偉

1西南交通大學機械工程學院 成都 610031 2西南交通大學牽引動力國家重點實驗室 成都 610031

0 引言

地鐵限界的計算方法經歷了不斷細化的發展過程。UIC505計算方法以系列標準的形式對車輛提出限界校核方法，但由于UIC標準提供的參考輪廓只適用于歐洲聯運鐵路，對于不是按照UIC標準修建的線路動態偏移量的設置方法是不同的，所以不能照搬其算法[1,2]。CJJ96—2003《地鐵限界標準》計算方法雖然考慮的因素較全面，但因無法準確考慮車輛一系二系振動變形量，且沒有車速與懸掛變形的對應關系，所以在確定計算參數時取值均趨向極限狀態可能會使計算結果偏于保守[3]。動力學計算方法沒有考慮到制造、安裝、維護誤差等因素，使得計算結果不完全切合實際[4,5]。在總結國內外限界計算方法研究基礎和經驗的基礎上，本文研究出一種基于Matlab編程與三維空間幾何變換理論[6]來獲取車輛動態包絡線的新方法，以某地鐵車輛為例，采用線路實測數據獲取車輛運行過程中的動態偏移量，再疊加CJJ96-2003《地鐵限界標準》[7]計算方法得到車輛的靜態偏移量，基于三維空間幾何變換理論推導方程，利用Matlab軟件編程來計算地鐵車輛在實際線路上運行的整個過程中車輛的動態輪廓，得到所有動態輪廓上控制點的坐標，并直接輸出車輛的動態包絡線，最終實現車輛動態包絡線的獲取。為了簡化計算，本文暫不考慮點頭、搖頭以及縱向位移對車輛動態輪廓的影響。

1 基于三維空間幾何變換的計算方程推導

圖1 位置示意圖

1.1 繞X軸旋轉變換

圖2 變化前后示意圖

1.2 三維平移變換

當空間向量發生平動時可以分解為在三個坐標軸方向的平動，假設在X、Y、Z方向平動值分別為TX、TY、TZ，則發生平動后其坐標只在每個方向的坐標值上分別疊加該方向的平動量。

1.3 三維空間幾何變換

當空間向量在三維空間內既平動又轉動時，假設先繞X軸轉動為再平動為其中轉角為φ，在X、Y、Z方向的平動量分別為。則變動前后的坐標關系為

且將矩陣稱為三維變換矩陣。

類似的在三維空間中存在某一剛體，其上有任意一點P= (x,y,z)，當剛體繞X軸旋轉φ角且發生平動后，該剛體上的P點變為P''點，兩點之間的關系也可以用上述公式表示。

2 計算截面與控制點

2.1 基準坐標系

軌道坐標的XOY平面位于軌面上，X軸在軌道中心線上并與車輛行進方向一致，Y軸與前者垂直并切于設計軌頂面，Z軸豎直向上，以X軸上任意一點為原點，構建笛卡爾右手直角坐標系O-XYZ[8]。

2.2 計算截面與車輛輪廓線坐標

為了描述整節車輛在線路上的運行行為，一般需在車輛上選擇多個截面輪廓來代表整個車輛，稱此截面為計算截面[9]。在計算車輛動態包絡線時，要先確定車輛的計算截面，車輛的計算截面應建立在車輛運動時車輛輪廓點發生最大位移的截面上，能體現出車輛的最大橫向、垂向動態偏移。

現選取車輛的4個計算截面,如圖3所示，分別為截面A-A（車體的端部截面）、截面B-B（過前轉向架中心的截面）、截面C-C（過車體中心的橫截面）以及截面D-D（過后轉向架中心的截面）。本文以某地鐵車輛為例按上述方法選取計算截面，考慮到篇幅，僅列出A-A斷面和B-B斷面上控制點坐標，如表1、表2所示。

圖3 截面位置示意圖

表1 靜止狀態下車輛A-A輪廓線控制點的坐標 mm

表2 靜止狀態下車輛B-B輪廓線控制點的坐標 mm

3 動態包絡線的獲取方法

3.1 獲取流程

本文以某地鐵車輛為例，研究獲取車輛的動態包絡線的新方法，獲取方法流程如圖4所示，通過線路實測數據獲取車輛運行過程中的動態偏移量再疊加由CJJ96-2003《地鐵限界標準》中的計算方法得到車輛的靜態偏移，基于三維空間幾何變換理論推導出方程，確定方程中的參數，利用Matlab軟件編寫函數來計算地鐵車輛在實際線路上運行的整個過程中車輛的動態輪廓及其上控制點的坐標，最后通過Matlab編寫程序直接輸出車輛的動態包絡線，最終實現車輛動態包絡線的獲取。

圖4 獲取方法流程圖

本文在研究車輛的動態包絡線問題時，先將計算截面簡化，如圖5所示。假設車體是絕對的剛體且質量分布均勻，采用前文所述的坐標系O-XYZ，要研究的問題可描述為已知剛性車體平衡狀態下的各點的坐標，車輛在線路上運行時產生振動，假設側滾角為φ，車輛發生的總的垂向偏移量為TZ，總的橫向偏移量為TY，總的縱向偏移量為TZ，車輛上存在一點Q，經車輛振動后Q(x，y，z)點最終變為點，兩者之間的關系可以用式（4）來表示。

圖5 計算截面簡圖

3.2 三維變換矩陣的參數確定

三維變換矩陣中的參數包括車輛在運行過程中的側滾角φ以及各方向上的偏移量TX、TY、TZ，本文暫不考慮車輛縱向位移帶來的影響所以TX的取值為零。

1）總的垂向偏移量TZ

假設計算截面A-A，B-B，N-N(任意截面)內由一系二系引起的車體的垂向偏移量分別為ZA、ZB、ZN，磨耗引起的垂向固定偏移量及各種安裝誤差等引起的垂向偏移量為Z、zi，其數值大小可參考CJJ96—2003《地鐵限界標準》得到。已知各輪對左側一系橫向坐標yaycyeyg大小為-965 mm，右側一系橫向坐標ybydyfyh大小為965 mm，前后轉向架左側二系橫向坐標ypyq大小為-925 mm，前后轉向架右側二系橫向坐標ymyn大小為925 mm，且已測得各輪對左右兩側一系引起的垂向偏移量za、zb、zc、zd、ze、zf、zg、zh，前后轉向架左右兩側二系引起的垂向偏移量zm、zn、zp、zq。經分析可得在任意截面N-N內，總的垂向偏移量為[4,10]

截面B-B內，由一系二系引起的車體的垂向偏移量zB，可分為兩部分，一是由一系引起的車體的垂向偏移量用zB1來表示，二是由二系引起的車體的垂向偏移量用zB2來表示，則有

分別求解zB1和zB2，因為車體是剛性的，所以對于截面B-B內任意點由一系、二系引起的車體的垂向偏移量zB1、zB2而言，其值與在一系、二系左右位置測量出的垂向偏移量有關，之間的關系可以近似為線性的

綜上可求出截面B-B內，由一系二系引起的車體的垂向偏移量ZB，令ZN=ZB再將其帶入式（5）中便得到截面B-B內，車體總的垂向偏移量TZ。用同樣的方法可以計算出截面D-D內總的垂向偏移量TZ。

已求出ZB、ZD后可計算除截面B-B、D-D外任意截面N-N內一系二系引起的車輛的垂向偏移量ZN，已知各截面相互平行如圖3所示，車輛定距a大小為12 600 mm，以截面C-C所在的位置作為X軸的原點，正方向與車輛行進方向一致，現在引入參數XN，表示任意截面N-N在X軸上的坐標，則有：

令XN=XA，XN=XC代入式（9）即可求得截面A-A，C-C內一系二系引起的車輛的垂向偏移量ZA、ZC。同樣在求其他任意計算截面內一系二系引起的車輛的垂向偏移量時，只需要將上式中的參數換成所要計算的截面的X軸坐標即可。綜上可求出任意截面N-N內由一系二系引起的車輛的垂向偏移量ZN，再將其代入式（5）中便得到任意截面N-N內總的垂向偏移量TZ。

2）車體總的橫向偏移量TY

車輛在線路上運行時，車輛相對于軌道發生的橫向偏移量包括：輪對相對于軌道的發生橫向偏移Y1，轉向架相對于輪對的橫向偏移Y2，車體相對于轉向架的橫向偏移Y3，以及磨耗引起的固定偏移量各種安裝誤差等引起的偏移量yi。設在任意截面N-N內輪對總的橫向偏移量為TY1，轉向架總的橫向偏移量為TY2，車體總的橫向偏移量為TY3，則有[4,10]

根據CJJ96—2003《地鐵限界標準》中規定的各個參數的值可以近似求得車輛總的橫向偏移量。

3）線路實測數據

運行過程中車體的側滾角以及由一系二系引起的垂向偏移是隨時間變化的。前者是在車體地板中心安裝角位移傳感器測得的，后者是在一系二系對應位置安裝位移傳感器測得的。為了使結果更精確更切合實際，在對實測數據進行采樣時，將采樣頻率設置為2 000 Hz，采樣時間500 s,最終得到106個采樣點，即za、zb、zc、zd、ze、zf、zg、zh、zm、zn、zp、zq，這些參數均取出個數值。由于實測數據較多,此處僅列舉抽取的側滾角φ，以及一位輪對左，一系引起的垂向偏移量za，如圖6、圖7所示。

圖6 側滾角φ

圖7 垂向偏移量za

4 基于Matlab獲取車輛動態包絡線的算例與分析

4.1 車輛動態包絡線的獲取

在利用推導出的式（4）來計算車輛偏移后各個計算截面內控制點的坐標時，需要很多輸入參數,如靜止狀態下計算截面輪廓線控制點的坐標x、y、z以及TX，TY，TZ，φ (其中x以及TX為零)，且本文抽取的線路實際測量的參數均有106個數值，將這些參數帶入前文所推導的公式中計算得到ZA、ZB、ZN也均有106個。由于輸入和輸出的參數數量龐大，本文利用Matlab軟件[10]自行編寫函數來計算車輛偏移后各個計算截面內控制點的坐標，輸出所有偏移的截面輪廓，并通過Matlab中的delaunayTriangulation和convexHull函數來實現車輛動態包絡線的獲取，計算流程如圖8所示。

圖8 計算流程圖

本文自行編寫的函數和腳本主要包括：

1）函數Bver，Dver，Nver，其作用是計算出各計算截面內由一系二系引起的車體的垂向偏移量；

2）函數Balt，Dalt，Nalt，其作用是計算出各截面發生偏移后各控制點的坐標；

3）Bplot，Dplot，Nplot，其作用是輸出各截面各個動態偏移輪廓；

4）Bmax，Dmax，Nmax，其作用是得到各截面各個控制點的極限坐標及最大的偏移輪廓，最終得到車輛的動態包絡線。

通過Matlab中的delaunayTriangulation和convexHull函數，繪制出這些離散點的最大輪廓，從而實現車輛動態包絡線的輸出[11]。此處僅列出程序中的部分代碼，具體內容如下所示：

DT=delaunayTriangulation(y,z);

%此函數基于列向量y和z中的點創建二維Delaunay三角剖分（y、z為控制點的坐標）

k=convexHull(DT)；

%此函數找到返回三角剖分的凸包的頂點

maxY=DT.Points(k,1)；

%各截面各個控制點的橫向的極限坐標

maxZ=DT.Points(k,2)；

%各截面各個控制點的垂向的極限坐標

plot(maxY,maxZ,‘r-’,‘LineWidth’,2)；

%畫出某計算截面的最大的偏移輪廓，即是所求的某計算截面內車輛動態包絡線

4.2 結果及分析

現將已得到車輛在某線路上運行的動態包絡線與CJJ96—2003《地鐵限界標準》中規定的限界進行對比校核。查找到A1型地鐵車輛在高架線曲線段設備限界坐標值輸入到Matlab中，最終得到利用本文所提出的方法獲取到的車輛動態包絡線與標準中的規定的車輛限界之間的關系，此處僅列舉計算截面A-A、B-B與標準對比的結果，如圖9、圖10所示。

圖9、圖10中黑色的線是車體靜態的截面輪廓，黑色的點是控制點的位置，綠色的線是計算截面各個偏移的輪廓，品紅色的線是利用本文所提出的方法獲取到的車輛動態包絡線，藍色的線是根據CJJ96—2003《地鐵限界標準》得到的車輛限界。從圖中可以看出本文所獲取到的車輛的動態包絡線遠小于標準中規定的車輛限界，即車輛滿足限界要求，可在線路上安全運行。

圖9 計算截面A-A的計算結果與對比圖

圖10 計算截面B-B的計算結果與對比圖

5 結論

本文基于三維空間幾何變換理論，采用線路實測數據并疊加靜態偏移量，利用Matlab軟件計算地鐵車輛在實際線路上運行的地鐵車輛動態輪廓，得到所有動態輪廓上控制點的坐標，并輸出車輛的動態包絡線，將其與CJJ96—2003《地鐵限界標準》中規定的限界進行對比分析，得到如下結論：

1)本文所提出的車輛動態包絡線獲取方法簡單便捷，可操作，易實現，采用線路實測數據并疊加靜態偏移使得最終結果更切合實際，可以更加有效地利用限界的有效空間，也為以后的車輛限界計算提供新思路。

2)利用Matlab自行編寫函數，能夠快速處理大量的數據，使得車輛動態包絡線獲取更迅速更便捷；通過編寫程序可直接輸出車輛動態包絡線，無需采用CAD繪制最大輪廓圖。

3)將本文所獲取的車輛動態包絡線與CJJ96—2003《地鐵限界標準》規定的限界進行對比分析，可知前者遠小于后者，即車輛滿足限界要求。