基于弧長參數化的曲面W-M分形插值

劉遠東， 王清輝， 劉 林， 熊 巍

（1. 華南理工大學設計學院，廣東 廣州 510006；2. 華南理工大學機械與汽車工程學院，廣東 廣州 510640）

有研究表明，機械加工表面的微觀粗糙形貌具有統計意義上的自仿射分形特性[1]，對機加表面可以像描述自然界的非規則結構一樣，使用分形幾何對表面微觀形貌進行表征建模。假設機加工表面微觀輪廓具有連續性，盡管這個假設在原子尺度時是不成立的，但對于機加工表面而言，一直到零長度尺度的極限情況，連續性假設都是成立的。據觀察，如果表面微觀輪廓被重復放大時，更加精細的結構會不斷出現，而且輪廓永遠是不光滑的，在任何點均不存在切線，所以輪廓是處處不可微性。另外，當輪廓被放大時，放大后的表面和原始表面的概率分布非常相似，因此大多數輪廓具有自仿射性。Mandelbrot于1982年在Weierstrass函數基礎上提出了一種分形曲線的函數表達形式，稱為Weierstrass- Mandelbrot函數（簡稱W-M分形）。分形幾何學中的W-M函數滿足連續性、處處不可微、自仿射性這些數學特征，因此Majumdar A用它來表征機械加工表面的微觀輪廓[2]。

Chen[3]等利用W-M分形模型和理論分析了微流道的層流表面，指出微流道的表面微觀形貌可以利用分形幾何的方法予以描述，而且自仿射的分形維數反映了表面粗糙形貌的微觀特征，自仿射分形維數越大，表面結構變化就越大，從而壓力損失也越大。Bobji[4]等利用用W-M分形理論建立了巖石表面結構的三維仿真模型，研究了分形維數D和均方根rms對巖石硬度的影響，從而利用W-M分形理論來進行巖石硬度的測量。李成貴[5]等指出工程表面的幾何特征與零件的摩擦、磨損、接觸、潤滑、磨削機理、金屬結構等機械或物理性能有著密切的聯系，通過計算W-M三維分形表面形貌的分形維數，并進行精密加工試件實驗，研究了各種分形維數與工程表面的形貌、加工參數、機械性能等的關系。He[6]等利用分形理論，對摩擦的合金表面進行取樣，通過分形參數D和G的分析，從而研究合金表面摩擦特性。呂建國[7]等研究了不同分形維數Ds的W-M分形曲面多重分形譜的變化規律，用分形參數定量表征分形曲面的分布特征，并指出多重分形譜很好地表征了W-M分形曲面的高度分布規律。張毅[8]將分形理論應用在掃描隧道顯微鏡測量的一組微鏡的表面形貌的分析中，指出W-M 函數可用于微鏡表面的建模和形貌分析。

然而，以上研究的都是基于直線上或平面上的W-M分形，在眾多的工程應用中，機械表面的微觀形貌往往是疊加在自由曲線或曲面上的。因此，對這些機械表面微觀形貌不能直接采用W-M分形進行模擬，而需要研究將W-M分形與自由曲線曲面合成的有效計算方法。本文針對具有自由曲面形狀的機械加工表面微觀形貌建模，研究將弧長參數化引入參數曲面上的W-M分形插值的計算方法。

1 W-M分形模擬

1.1 W-M二維曲線模擬

分形幾何學中W-M函數的數學表達式為[9]

式中，x表示輪廓的位置坐標，D為分形維數，它描述函數z(x)在所有尺度上的不規則性，但不能確定z(x)的具體尺寸；G是特征尺度系數，反映 z(x)幅值大小，它決定 z(x)的具體尺寸；γn表示輪廓的空間頻率，γ為大于1的常數，對于服從正態分布的隨機輪廓，γ=1.5可適用于高頻譜密度及相位的隨機性，由于粗糙度輪廓是非穩定的隨機過程，輪廓結構的最低頻率與粗糙度樣本長度關系為 γn1≈1/L1；n1是W-M函數的初始項，為整數；L1為粗糙度樣本取樣長度。

圖1是基于上述W-M分形公式得到的一段分形曲線輪廓。

圖1 W-M二維分形插值曲線

1.2 W-M三維曲面模擬

W-M 函數實際上是空間不同層次不同頻率的正弦/余弦曲線疊加而成，其中一種常用的產生隨機分形表面的典型函數由式（2）表示[10]

其中，C是表面的高度系數，Ds是空間分形維數且有 Ds=D+1，D是截面的分形維數，Φ是一系列0到π之間的隨機相位，γ是控制表面起伏密度的參量，M是空間起伏的數量，L為擬圖像的大小，G則是調整粗糙度的一個系數。對于一個完全分形的表面，有nmax→∞ 。圖2是在Object ARX環境下運用W-M三維模擬公式得到的一張分形插值曲面模型。

圖2 W-M三維分形插值曲面

2 弧長參數化

2.1 弧長參數化的定義

弧長參數化是指基于曲線自身的累積弧長建立參數方程。在許多實際問題中，需要對曲線進行弧長參數化。例如，數控加工中要求實現刀具運動軌跡的弧長到空間位置的“時空”轉換，這里要求弧長作參數，這樣只要給出弧長即可立刻計算刀具的位置。

以 sin(t)曲線的弧長參數化為例，插值時以均勻遞增的累積弧長作為W-M分形插值的參數變量，可得到W-M分形曲線與sin(t)曲線疊加的分形模型。圖3是利用弧長參數化得到的 sin(t)曲線上的W-M分形曲線例子。

圖3 sin(t)曲線上的W-M分形曲線

2.2 弧長參數化必要性

對參數方程為P(t)的曲線進行W-M分形插值時，如果直接以均勻遞增的參數 t作為 W-M分形插值的參數變量，必然導致W-M分形插值沿曲線表面分布的不均勻。在曲線參數密度大的地方，分形插值點會相對密集；在參數密度小的地方，W-M分形插值點相對稀疏。

圖4(a)是直接基于曲線參數變量驅動得到的一條NURBS曲線上的W-M分形插值曲線，該曲線左邊參數密度大處，W-M分形插值點密集；而曲線右邊參數密度較小，造成W-M分形插值分布稀疏。圖4(b)是采用弧長參數作為驅動變量得到的NURBS曲線上的W-M分形插值曲線，插值點是均勻的，它不隨曲線的參數密度變化而變化。由于機械表面的微尺度形貌特征是獨立與其宏觀幾何形狀而存在的，而更多地受到表面微分幾何特征的影響，因此，采用弧長參數作為驅動變量對于實現參數曲線上的W-M分形插值是必要的。

圖4 兩種參數驅動生成的W-M分形插值曲線對比

3 曲面W-M分形模擬

3.1 弧長參數化的實現

本文以AutoCAD環境下的一張NURBS拉伸曲面為例，說明弧長參數化的實現方法。假設，W-M分形插值模擬分辨率為m_nCol×m_nRow個插值點，且曲線的拉伸方向插值點個數為m_nCol，弧長參數化方法如下：

步驟1將曲線的參數區間(ts，te)分為n等份（n=m_nRow×Pric，Pric為設定的弧長參數化精度），利用微分法求出每一ti(i=0,1,…,n-1)所對應累積弧長si，故曲線總弧長為sn-1。

步驟2反求出m_nRow個等弧長點ks～對應的參數tk，其中m_nRow-1。讓ks ～在弧長區間 ( s0,sn-1)內由s0開始搜索，若判得，則認為si≈ks ～，保存此時的參數值tk。由于本插值問題不同于一般的弧長參數化，對每一個進行搜索時，不必均從s0開始搜索，因為弧長S(t)是單調遞增函數，所以，若上一次迭代得(k-1)對應的值為sj，那么ks ～在區間(sj,sn)從sj開始搜索，這將大大減少弧長參數化的算法實現時間。

步驟3前面已求得等遞增的累積弧長對應的參數值tk，由Object ARX的求值函數可求出tk對應的坐標Pk和曲線法向量Vk，再求出W-M分形插值的高度值Hk，根據Pk、Vk和Hk可求出最后的W-M分形插值點。對于二維W-M分形插值，把ks ～作為x值代入式(1)求Hk；對于三維W-M分形插值，若拉伸曲面沿y方向拉伸，則y方向的值不變，ks ～代替x，代入式(2)求Hk。

3.2 W-M分形模擬對比

對于曲線上的二維 W-M 分形插值模擬對比，圖4給出了模擬結果，其中插值點個數為40。對于曲面上的三維 W-M 分形插值，圖 5(a)是一張NURBS拉伸曲面，圖5(b)顯示了直接以曲面參數進行插值得到的分形曲面，而圖5(c)則是以弧長參數化方式得到的分形插值曲面，其中插值分辨率均為45×15。從圖形對比可以觀察到，采用弧長參數驅動生成的分形插值曲面其分形形貌分布均勻。

圖5 兩種參數驅動生成的W-M分形插值曲面對比

4 兩向異性的W-M分形

在對上述表面微觀形貌進行模擬時，我們假設其微觀形貌特征的分布是各向同性的。但事實上，很多機械加工過程產生的表面其微觀形貌特征呈很強的方向性[11]。工程一般認為，各向異性表面的加工紋理，就是加工痕跡的線條方向。紋理的種類有平行、垂直、相交、微粒狀、圓形、放射狀[12]。這里我們先以沿x方向和y方向的兩向異性情況為例研究其數學模型的表示。

本文提出以x方向的二維W-M分形和y方向的二維W-M分形相疊加的方法來實現這種兩向異性的W-M分形模擬，公式如下

其中，Gx和Dx為沿x方向的分形比例系數和分形維數；Gy和Dy為y方向的分形比例系數和分形維數。

由式(3)可知，當x或y為一定值時，z(x,y)演變為二維W-M分形輪廓線。這也表明，利用該式模擬的曲面，在x方向或y方向上的任意輪廓線均具有分形特征。

圖6(a,b)分別是由式(3)進行W-M插值得到曲面，兩插值曲面的尺度范圍不同，其中圖6(a)為 100×100，圖 6(b)為 10×10。從圖中可看出，曲面的x方向和y方向是異性的，且兩曲面具有自相似性。圖 6(c)表明當 Gx和 Gy相等、Dx和Dy相等時，插值曲面在x和y兩方向同性。圖6(d,e)則顯示在 Gx、Gy相等的情況下，由于 Dx和 Dy不同得到x方向和y方向兩向異性的分形曲面。

圖6 兩向異性W-M分形模擬曲面

通過前面提出的弧長參數化方法以及各向異性W-M分形插值公式，可以在參數曲面上合成具有的各向異性W-M分形插值曲面。圖7展示了兩張在NURBS拉伸曲面上合成的各向異性W-M分形曲面。

圖7 曲面的兩向異性W-M分形模擬曲面

5 結 論

本文論述了在參數曲面上合成W-M分形在機械加工表面微觀形貌建模中的意義，指出了弧長參數化在該建模過程中的必要性，并提出了具體實現方法。通過該方法，利用 AutoCAD ObjectARX開發包實現了曲線和曲面的W-M分形插值模擬。在此基礎上,提出了一種兩向異性W-M 分形曲面的數學建模方法，實現了在參數曲面生成兩向異性W-M分形插值模擬。

[1] 董 霖,王 艷,王 霜. 機械加工表面分形特征的數學表征[J]. 機械,2001,28: 52-53.

[2] Majumdar A,Tien C. Fractal characterization and simulation of rough surfaces [J]. Wear,1990,136:313-327.

[3] Chen Y P,Zhang C B,et al. Role of surface roughness characterized by fractal geometry on laminar flow in microchannels [J]. Physical Review,2009,80(026301):1-7.

[4] Bobji M S,Shivakumar K,et al. Influence of surface roughness on the scatter in hardness measurements-a numerical study [J]. International Journal of Rock Mechanics and Mining Sciences,1999,36: 399-404.

[5] 李成貴,董 申. 三維表面形貌的分形維數計算方法[J]. 航空精密制造技術,2000,(8):36-40.

[6] He L,Zhu J. The fractal character of processed metal surfaces [J]. Wear,1996,208(197): 17-14.

[7] 呂建國,戴結林,宋學萍. Weierstrass-Mandelbrot分形曲面的多重分形譜[J]. 功能材料,2008,(9):1574-1576.

[8] 張 毅,羅 元,徐曉東. 基于分形理論的 MEMS光開關微鏡分析與模擬[J]. 半導體光電,2006,(6):672-674.

[9] Majumdar A,Bhushan B. Fratal model of elasticplastic contact between rough surfaces [J]. Journal of Tribology,1991,113: 1-13.

[10] Yan W,Komvopoulos K. Contact analysis of elastic-plastic fractal surfaces [J]. Journal of Applied Physics,1998,84(7): 3617-3624.

[11] 喻志清,戴澤西,陳志祥. 外圓磨削表面輪廓的分形行為研究[J]. 華中理工大學學報,1998,(3):49-50.

[12] Thomas T R,Rosen B G,Amini N. Fractal characterization of the anisotropy of rough surfaces [J]. Wear,1999,232: 41-50.