無奇異性邊界積分法在非平滑物體之應用

1

行政院國家科學委員會專題研究計畫成果報告

無奇異性邊界積分法在非平滑物體之應用

Nonsingular boundar y integr al methods for nonsmooth bodies

計畫編號：NSC 89-2611-E-002-046

執行期限：89 年 8 月 1 日至 90 年 7 月 31 日

主持人：黃維信 國立台灣大學造船及海洋工程學系

計畫參與人員:洪立萍 國立台灣大學造船及海洋工程學系

計畫參與人員:葛家豪 國立台灣大學土木工程學系

一、中文摘要 本文旨在建立無奇異性邊界積分式的 推導。利用高斯通量定理以及將積分式中 之 奇 異 性 核 函 數 加 入 一 項 使 其 奇 異 性 消 除，再將加入函數的解析積分值減回原方 程式。經過去除奇異性的步驟後，無奇異 性邊界積分方程式即可直接在積分邊界離 散。離散方式為直接將分布於積分邊界的 積分點視為節點，其滿足邊界條件並組成 影響係數矩陣，以求的流場中之速度勢。 本文中並探討邊界積分法在平滑邊界流場 與非平滑邊界流場中之特性。 關鍵詞：邊界積分法、勢流理論、非奇異積 分方程 Abstr act

This study presents a general non-singular boundary integral equation for potential flow. The singularities of Green’s function and its normal derivative are desingularized theoretically, using a subtracting and add back technique, and Gauss flux theorem, so that any numerical integration methods can be applied directly evaluate all the integrals without any difficulty. The collocation points are chosen to be the nodes of numerical integration. Several numerical examples are calculated

with smooth boundary and non-smooth boundary by boundary integral method.

Keywor ds: boundary integral equation;

potential flows; non-singularity integral equation 二、緣由與目的 勢流場的相關問題常常出現在工程與 物 理 問 題 上 ， 而 勢 流 場 的 控 制 方 程 式 為 Laplace 方程式。以往求解微分方程所運用 數值方法，如有限差分法(Finite Difference Method)、有限元素法(Finite Element Method) 等，都必須將計算域做空間的離散。本研究 使用另一種解法，不需要做空間的離散，即 是將 Laplace 方程式利用格林定理轉換為邊 界積分方程式來進行求解。轉換為邊界積分 方程形態後，其計算域即降低一階，也就是 若為三維的問題，只須對其邊界曲面進行運 算，而若為二維的問題，只須對其邊界曲線 作運算。這個特性在數值計算上不但減低了 相當大的計算量及記憶體空間，也縮短了求 解的時間。 傳統的邊界元素法(Boundary Element Method)常運用於計算邊界積分方程式。其做 法為先在邊界上進行分割，形成邊界元素的 組合，並引入形狀函數(shape function)為元 素中之每個物理量的變化基礎，再針對每一

2 個元素進行積分的離散。在此，形狀函數的 引入使得計算積分的邊界幾何形狀失真，若 原本之外形有轉折或曲率變化較大者，則需 要更多的元素來近似其形狀以減少空間上的 離散誤差。再者，形狀函數對於未知的函數 限制其變化的趨勢，亦使原函數喪失部份原 有的特性。 本 研 究 計 劃 以 邊 界 積 分 法 (Boundary Integral Method) 直接進行邊界積分方程的 數值計算。邊界積分法主要是先將邊界積分 方程式以代數的技巧消除其奇異性，而得到 一無奇異性之邊界積分方程式。無奇異性的 邊界積分法即可直接在邊界上分布積分點的 方式，離散積分式以求得積分值。在計算的 過程中並不存在任何的形狀函數，對於未知 的函數亦無任何限制，所以邊界積分法改善 了 BEM 中引入形狀函數所可能產生的誤 差，並且在計算的過程中也因少了形狀函數 的運算而減少了一部份的計算步驟。因此在 使用邊界積分法計算邊界積分式時，只要能 掌握離散求解積分值時所產生的誤差就可以 得到快速又精準的計算結果。 由於引入了格林函數，邊界積分方程式 中存在著因格林函數及其法方向微分量所產 生 之 有 奇 異 性 的 積 分 核 函 數 (singular kernel)。在本研究中，我們將以 subtracting function[1-6]的概念消除積分式的奇異性。在 平滑邊界的問題中，我們應用高斯通量定理 (Gauss flux theorem)消除格林函數之法方向 微分量的奇異性。對於格林函數的奇異性， 我們應用 Saranen [7]的方法，將積分邊界映 射至一個單位圓上。在靠近奇異點時，利用 L’Hopital’s Rule 即可計算出其極限值。在非 平滑邊界的問題中，我們先將邊界分割為數 段平滑的曲線。相同地，我們應用高斯通量 定理消除格林函數之法方向微分量的奇異 性。再將含有奇異點的邊界映射至一直線段 上，並利用 L’Hopital’s Rule 計算出在靠近奇 異點時之極限值。經過以上消除奇異性的步 驟，原來的邊界積分方程式即轉換為無奇異 性邊界積分方程式並可以直接利用數值積分 方法將積分邊界離散，以求得流場之速度勢 及其方向微分量的代數方程式進而求解其 值。 三、結果與討論 3.1 Smooth boundar y

The first test case is a Dirichlet problem with a boundary condition, φ = 1, on an ellipse, x2+16y2=1. When the trapezoidal rule is chosen as the quadratural rule, the nodes are uniformly distributed along the polar angle. When the Gaussian quadrature rules are applied, the ellipse is split into upper and lower parts. The integration points are distributed along the polar angle in the range [0, π] on each semi-ellipse in order to maintain the property of symmetry. Though both methods with double precision almost give exact solutions as shown in Fig. 1, the trapezoidal rules relatively give better results. It is well known that the trapezoidal rules have very fast convergent speed for periodic functions in the numerical quadrature. In this case, since the accuracy is up to the limit of computer, the primary error is due to the error accumulation of computing. Therefore, the error increases when the number of computing nodes increases.

3.2 Nonsmooth boundar y

A unit square is chosen as the test case for the nonsmooth boundary shape. In this geometry, there are four discontinuous normal vectors at corners. So, the normal derivative of potentials and the normal derivative of Green's function are also discontinuous at corners. The solution of test cases is chosen as before φ=sin(3x) cosh(3y). In these cases, we avoid trapezoidal rules for the sake of the discontinuity of functions at corners. The first one is the Dirichlet problem and the results are shown in Fig. 4. The results indicate that the convergent

3

speed for both the constant element method and the present method with Gaussian nodes is of the same order O(N-1/2), when the number of nodes is large enough. For accuracy, the Gaussian formulas are better than the constant element method in the root-mean-square sense. However, Fig. 4 cannot show the detailed variation of errors. Therefore, a typical case with sixteen Gaussian points or sixteen constant elements along each side of the square is presented in Fig. 3. As expected, it shows that most errors are concentrated near corners, which is known as corner singularities or geometrical singularities. Since the errors near corners are much bigger than other places, they dominate the magnitude of root-mean-square errors. It makes the convergence speed become slow in the least-square sense, when the nodes increase as shown in Fig. 2. From Fig. 3, it also confirms that the present method is much better than the constant element method for the error at the nearby location. However, the difference of root-mean-square errors between both methods cannot present this. It is because the Gaussian formula contains more nodes near corners, while the constant element method has few. Besides that, the end points of Gaussian formulas are closer to corners. Therefore, the difference of root-mean-square errors between these two methods is diminished after avera-ging by the number of nodes.

四、參考文獻

1. Kantorovich LV, Krylov VI. Approximate methods for higher analysis. Interscience, 1958.

2. Landweber L. and Macagno M. ,‘Irrotational flow about ship forms ,’ IHHR report , Iowa , No.123 , 1969.

3. Hwang W.S. ,‘Hypersingular boundary integral equations for exterior acoustic problems,’ J. Acoust. Soc. Am. 1997 ; 101(6) : 3336-3342.

4. Hwang W.S. and Huang Y.Y. ,‘Non-singular direct formulation of boundary integral equations for potential flows,’ Int. J. Numer. Mech. Fluids 1998 ; 26 :627-635.

5. Chang J.M. ,‘Numerical studies on desingularized

Cauchy’s formula with applications to interior potential problems ,’Int. J. Numer. Mech. Engng.1999; 46: 805-824.

6. Yang S.A. , ‘On the singularities of Green’s formula and its normal derivative with an application to surface-wave-body interaction problems ,’ Int. J. Numer. Mech. Engng. 2000 ; 47(11) : 1841-1864.

7. Saranen J. , ‘The modified quadrature method for logarithmic-kernel integral equations on closed curves ,’ Journal of Integral Equations and Applications 1991 ; 3(4) : 575-600.

五、成果自評

本研究內容與原計畫相符、並達成預期 目標、研究成果兼具學術及應用價值、適合 在學術期刊發表。

4 2 3 5 2 3 10 100 N 1E-15 1E-14 1E-13 1E-12 E 2

Fig. 1 Root-mean-square errors via number of total nodes on an ellipse (a:b=1:1/4) with the boundary condition f=1; ì (Gaussian quadrature), â (trapezoidal rule) 6 8 2 4 6 8 2 4 6 8 10 100 N 1E-1 1E+0 1E+1 E 2

Fig. 2 Root-mean-square errors via number of total nodes on a unit square with the Dirichlet boundary condition φ=sin3x cosh3y; è (constant element), ì (Gaussian quadrature)

Fig. 3 The distribution of absolute errors on each side of a unit square with the Dirichlet boundary condition φ=sin3x cosh3y; è (constant element), ì (Gaussian quadrature); distribution on (a) y=0, (b) x=1, (c) y=1, (d) x=0 0.0 0.5 1.0 x 1E-4 1E-3 1E-2 1E-1 1E+0 0.0 0.5 1.0 y 1E-4 1E-3 1E-2 1E-1 1E+0 1E+1 0.0 0.5 1.0 x 1E-4 1E-3 1E-2 1E-1 1E+0 1E+1 0.0 0.5 1.0 y 1E-4 1E-3 1E-2 1E-1 1E+0 1E+1 (a) (b) (c) (d)