參考資料

G.R.Liu Y.T.GU著 王建明 周學軍譯 《無網格法理論及程序設計》

數值實現

Matlab 2019a

前情回顧

形式主義的居士：無網格法理論與Matlab程序設計（1）——概述?zhuanlan.zhihu.com地球物理局 地震波場模擬實驗室 無網格組 地球物理局 基建處 數值計算科聲明：# 系列寫作內容首先符合本人的研究需要，不會優先照顧讀者體驗。 # 僅供學習和參考，禁止轉載。

---

傳統RPIM形函數

為避免采用多項式基PIM所引起的奇異性問題，徑向基函數（RBF）被采用以形成徑向基點插值法（RPIM）形函數用于無網格弱式法（GR Liu和Gu，2001c；Wang和Liu，2000；2002a，c）。RPIM形函數在無網格弱式法和強式法中均得到廣泛應用。

添加多項式的RPIM插值可表示為：

式中，

為徑向基函數（RBF）， 為徑向基函數的個數， 為空間坐標 的單項式， 為多項式基函數的個數。如 ，則為單純的徑向基函數；否則為添加了基函數的徑向基函數。 和 為待定常數。

在徑向基函數

中，僅有表示計算點 與節點 之間距離的變量

常用的徑向基函數（RBF）有幾種，這些RBFs的性質均已得到深入研究。最常用的4種分別為復合2次（MQ）、Gaussian（Exp）函數、薄板樣條（TPS）函數以及對數型徑向基函數。如下所示。在使用這些RBFs時，要確定幾個形狀參數，通?？筛鶕o定問題的數值實驗來確定。

復合2次（MQ）：

Gaussian（EXP）：

薄板樣條（TPS）：

對數型：

注：

為一與計算點 的局部支持域中的節點間距有關的特征長度；該長度通常等于局部支持域中的節點平均間距。

另外還存在一類所謂緊支徑向基函數（CSRBFs）。傳統RBFs對于

小于或等于某一定值，為嚴格正定的，故可形成理想的 階光滑。在CSRBF中，用一個形狀參數 決定該CSRBF的局部支持域尺寸。當 時，其值視為零。然而CSRBFs較傳統RBFs在表面擬合和求解力學問題時沒有明顯優勢。

Wu-C2：

Wu-C4：

Wendland-C2：

Wendland-C4：

Wendland-C6：

為局部支持域尺寸。

（1）中的多項式項并不總是必需的。但采用單純的徑向基函數的RPIM形函數不能通過標準分片試驗。另外添加多項式可以提高解的精度、使受形狀參數影響的敏感性降低、對于某些徑向基函數可提升穩定性。

為了確定（1）中的

和 ，需要形成計算點 的支持域，其中包括 個場節點。使（1）滿足計算點 周圍的 個節點值以確定系數 和 。這將產生 個線性方程，一個節點對應一個方程。這些方程可以表示為矩陣形式：

式中的函數值向量

為

徑向基函數的力矩矩陣

為

多項式力矩矩陣

為

徑向基函數的系數向量為

多項式系數向量為

在（5）中，

的 定義為

然而（3）中有

個變量，可以使用下面 個約束條件添加 個方程

聯立（12）和（19），得到如下矩陣方程

式中

因為矩陣

對稱，故矩陣 也是對稱的。求解（11）可得到

可將（1）重寫為

利用（23）可得到

RPIM形函數可以表示為

對應于節點位移向量的RPIM形函數

可表示為

（25）可重寫為

可方便地得到

的導數

其中

表示坐標 或 ，逗號表示對其后的空間坐標求偏導數。

須指出，對于任意的節點分布

通常是存在的。另外由于（1）中所使用的多項式階數相對較低，故通常情況下，RPIM局部支持域中使用少量節點（對于二維問題通常是10-40）不會引起奇異性問題。

隨著節點數增加，力矩矩陣

的條件數將會變差?？赏ㄟ^無網格全局配點法，即使用整個問題域中的所有節點構造計算公式，觀察到上述結論。全局配點法的一個特性是有可能構造對稱公式，在此不討論。

---

優缺點

使用局部支持域和RBFs基構造PIM形函數有以下優點：

• 使用RBFs可有效克服多項式PIM的奇異性問題。
• RPIM形函數是穩定的，故可適應任意不規則的節點分布。
• 便于構造三維RPIM形函數，因為在RBF中僅包含距離 變量。對于三維插值可方便地將距離表示為

• RPIM形函數較MLS形函數能更好地適應流體動力學問題。

但RPIM也存在一些缺點，如：

• RPIM形函數在解決固體問題時通常較MLS形函數和PIM形函數求解精度低。
• 一些形狀函數需仔細決定，因為它們將影響在無網格法中使用RPIM形函數的求解精度及功效。
• RPIM形函數通常較多項式PIM的計算費用高，因為其要求局部支持域中包含更多的節點。

---

形函數性質

RPIM形函數性質為：

（1）

函數性

RPIM形函數具有Kronecker

函數性質。

（2）單位分解性

RPIM形函數擁有單位分解性質，即

如果將線性多項式添加到基中，則在基函數中含有常數項。該單位分解性質可利用RPIM形函數所擁有的的再生性質方便地得到證明。

對于CSRBF，如果采用純粹的RBF，該單位分解性可很容易地得到證明，因為在CSRBFs中明顯包含常數項。然而在某些RBFs中沒有明顯的常數項，如MQ-RBF。對于這些RBFs需采取必要措施使其明顯包含常數項。

可將一個具有各階連續導數的函數展開為無窮Taylor級數。如對于MQ-RBF，在

附近的Taylor級數展開為

由式（32）可以清楚地看到其中包含一個常數項，因為MQ-RBF中的

。該常數基的存在可使RPIM形函數通過其再生性質再生出常數場。注意在局部域中再生一常數場的條件是在RPIM中利用各種RBFs獲取精確解所必需的，因為式（3.92）可展開至無窮階項。故式（31）在數值驗證時只能近似而不能精確滿足，由于數值舍棄誤差的存在使得RBF的計算總是對應于有限項的Taylor級數展開。

注意，TPS-RBF和對數型RBF不滿足

的條件，故在使用TPS-RBF和對數型RBF時須添加多項式項以確保其RPIM形函數擁有單位分解性。

（3）緊支性

RPIM形函數是緊支的，因為它是通過一緊支域的節點構造的，在該支持域以外不使用該形函數或認為其值為零。

（4）連續性

RPIM形函數通常擁有高階連續性，因為其徑向基函數是高階連續的。

（5）再生性

至少含多項式項的RPIM能確保精確再生線性多項式。

注意，有些RBFs不具有線性再生性，即其RPIM在沒有添加線性多項式項時不能再生任一線性場函數。例如在MQ-RBF中，由式（32）所示的其Taylor展開中不含線性項，因為

。這可能是使用MQ-RBFs所得到的 收斂性不佳的主要原因，故有時需在RPIM中添加線性多項式項以改善其性能。

（6）協調性

在使用RPIM形函數時，如采用局部支持域，則總體域上的協調性無法保證。當節點進入或離開移動支持域時，其場函數的近似式可能是不連續的。

---

數值實現

確定支持域

對于一個二維域

，計算點的支持域常用圓形和矩形。矩形便于構造和應用，在此我們使用矩形域。

矩形支持域的尺寸分別由

方向的 和 決定，即

其中

和 分別為該支持域沿 和 方向的無量綱尺寸。為簡化起見，常使用 。 和 為插值點 附近沿 和 方向的節點間距。如節點均勻分布，則 可簡化為沿 方向的兩相鄰節點之間的間距，而 可簡化為沿 方向的兩相鄰節點之間的間距。如節點非均勻分布， 和 可由求支持域中平均節點間距確定。

2. 徑向基函數形狀參數

包括復合2次(MQ)-RBF、Gaussian(EXP)-RBF以及薄板樣條(TPS)-RBF。RBFs中的參數需要事先確定。

（1）對于MQ-RBF，有兩個形狀參數

和 。在傳統RBF中常取 。對于固體力學和流體力學， 或 時結果較好。節點間距 可表示為：

式中

和 分別為定義在（1）中的沿 和 方向的節點間距。

（2）對于EXP-RBF，僅含一個形狀參數

，該值通常小于 。

（3）對于TPS-RBF，僅含一個形狀參數

。

形狀參數將對RBFs的性能產生影響，通常沒有最佳理論值。

3. RPIM形函數計算

RPIM形函數可表示為：

RPIM形函數在解決固體問題時通常較MLS形函數和PIM形函數求解精度低。

采用線性方程求解器從而避免直接對

求逆。（3）可重新寫為：

故可得到

或

用標準線性方程求解器求解（6），可以直接得到RPIM形函數而不需計算

。而且利用（6）可以計算RPIM形函數的導數

或

2階導數為

故利用標準線性方程求解器分別求解（8）和（9）可獲得RPIM形函數的導數。

4. 子程序流程圖

---

RPIM形函數算例

通過一個支持域中包含25個節點的算例，說明傳統RPIM形函數及其導數的性質。

該

個節點規則而均勻地分布在一矩形域： 和 。節點分布如上圖和表所示。分析采用 個采樣點用于計算形函數和繪圖。

我們給出一簡單的主程序，采用影響域替代支持域。調整不同場節點影響域大小，以使所有

個場節點均包含在每個計算點的插值計算中。

RPIM_main.m（程序疑有誤，建議對照原書檢查，希望得到指正）

clear

RPIM_ShapeFunc_2D.m（程序疑有誤，建議對照原書檢查，希望得到指正）

%%%%%% 參考《無網格法理論及程序設計》%%%%%%%%%%%%%

Compute_RadialBasis.m（程序疑有誤，建議對照原書檢查，希望得到指正）

% 參考《無網格法理論及程序設計》編制，源程序為Fortran語言

GaussEqSolver_Sym.m（程序疑有誤，建議對照原書檢查，希望得到指正）

% 參考《無網格法理論及程序設計》編制，源程序為Fortran語言

---

輸出結果（與原書結果并不匹配）

Node x y Phi dPhidx dPhidy dPhidxx dPhidxy dPhidyy 1 -1.000 -1.000 2.09762 0.71158 0.83018 0.29273 -0.35030 0.184302 -1.000 -0.500 1.80278 0.75895 0.56921 0.26816 -0.27322 0.427543 -1.000 0.000 1.61245 0.79736 0.26579 0.24201 -0.14082 0.617534 -1.000 0.500 1.56525 0.80828 -0.06736 0.23352 0.03667 0.670515 -1.000 1.000 1.67332 0.78419 -0.39209 0.25162 0.20093 0.553026 -0.500 -1.000 1.85742 0.43705 0.87410 0.50510 -0.23852 0.147317 -0.500 -0.500 1.51658 0.47849 0.61520 0.52706 -0.20122 0.424858 -0.500 0.000 1.28452 0.52085 0.29763 0.54222 -0.11535 0.678169 -0.500 0.500 1.22474 0.53576 -0.07654 0.54568 0.03138 0.7608810 -0.500 1.000 1.36015 0.50492 -0.43279 0.53742 0.15763 0.5862111 0.000 -1.000 1.73205 0.12872 0.90103 0.63293 -0.07464 0.1210912 0.000 -0.500 1.36015 0.14426 0.64919 0.70631 -0.06756 0.4173213 0.000 0.000 1.09545 0.16625 0.33250 0.80828 -0.04595 0.7393514 0.000 0.500 1.02470 0.18080 -0.09040 0.87447 0.01478 0.8966315 0.000 1.000 1.18322 0.15653 -0.46960 0.76349 0.05753 0.6100716 0.500 -1.000 1.74642 -0.19238 0.89777 0.61753 0.11075 0.1244117 0.500 -0.500 1.37840 -0.21491 0.64472 0.68327 0.09925 0.4185918 0.500 0.000 1.11803 -0.24495 0.32660 0.76751 0.06532 0.7294019 0.500 0.500 1.04881 -0.26149 -0.08716 0.81203 -0.01987 0.8650120 0.500 1.000 1.20416 -0.23209 -0.46418 0.73196 -0.08334 0.6069421 1.000 -1.000 1.89737 -0.49496 0.86617 0.46713 0.26524 0.1545222 1.000 -0.500 1.56525 -0.53885 0.60621 0.47799 0.22002 0.4260423 1.000 0.000 1.34164 -0.58123 0.29062 0.48109 0.12273 0.6651824 1.000 0.500 1.28452 -0.59526 -0.07441 0.48042 -0.03296 0.7399525 1.000 1.000 1.41421 -0.56569 -0.42426 0.48083 -0.16971 0.57983

---

總結

以上是生活随笔為你收集整理的c2064 项不会计算为接受0个参数的函数_无网格法理论与Matlab程序设计（6）——传统径向基点插值（RPIM）形函数...的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。