流体解析ソフトウェア Particleworks
  • Home
  • Home
  • Case Examples
    • Particleworks解析事例
    • Particleworks Case Examples
  • Technical
    • 粒子法・MPS法
    • 技術コラム >
      • DX時代の製品開発プロセスとCAEの重要性 >
        • 第1回 序 略歴とコラム紹介
        • 第2回 DXとデジタルエンジニアリング
        • 第3回 製品開発プロセスの目指す姿
        • 第4回 DX時代のCAE
        • 第5回 評価CAEの概要と課題
        • 第6回 評価CAEの課題解決手法
        • 第7回 企画CAEの概要と課題
        • 第8回 企画CAEの運用と応用
        • 第9回 設計CAEの概要と課題
        • 第10回 設計CAEの課題解決の進め方
        • 第11回 開発プロセス運用の仕組み作り
        • 第12回 まとめと変革の時代に求められるエンジニア像
      • 粒子法のいま、そして未来へ >
        • 第1回 粒子法のいま
        • 第2回 SPH法におけるカーネル近似とカーネル関数の条件
        • 第3回 SPH法における空間離散化
      • 粒子法の非圧縮条件とは
      • 粒子法入門 >
        • 第1回 粒子法って何?
        • 第2回 粒子法は、他の方法とどう違うか
        • 第3回 粒子法の大きさと質量について
        • ​第4回 「粒子の動かし方」と「加速度の求め方」について
        • ​第5回 計算時間を短縮する方法について
    • Technical Column >
      • Growing the particle method, and its present state >
        • 1. Present State of the Particle Method
        • 2. Kernel Approximation and Kernel Function Conditions in the SPH Method (Preparation for Spatial Discretization)
      • Incompressibility of the particle method
      • Introduction to the particle method >
        • 1. What is a particle method?
        • 2. In what ways is the particle method different from other methods?
        • 3. Mass and volume of particles
        • 4. How to move particles and how to calculate accelerations of particles
        • 5. How to shorten the simulation time
    • 粒子法用語集
    • Particle Method Glossary
    • 参考文献・ウェブサイト
    • Reference Book/URL
    • 論文・講演
  • Contact
    • 導入の流れとライセンス形態
    • Particleworks / GranuleworksプリインストールGPU搭載ワークステーション
    • 開発元・パートナー
    • Developers, Partners
    • お問い合わせ
    • Contact Us

2. Kernel Approximation and Kernel Function Conditions in the SPH Method (Preparation for Spatial Discretization)

This section explains a spatial discretization scheme in the Smoothed Particle Hydrodynamics (SPH) method. In the SPH method, a continuous physical quantity is approximated by a kernel function, and then the kernel approximation becomes a discretize form with particle quantities in the vicinity of the target point. For example, a scalar function \( \phi\) of the field defined on a three-dimensional (or two-dimensional) region \( \Omega\) is expressed in the form of a volume integral using the Dirac delta function \( \delta\) at an arbitrary point \( \boldsymbol{x}\) in the region.

$$ \phi(\boldsymbol{x})=\int_{\Omega} \phi(\boldsymbol{\xi}) \delta(\boldsymbol{\xi}-\boldsymbol{x}) \mathrm{d} \boldsymbol{\xi} $$

Here, the symbol \(\int_{\Omega}(\cdot) \mathrm{d} \boldsymbol{\xi}\), which denotes the volume integral, refers to the integration of the volume, where the integrant is a function of the variable \(\boldsymbol{\xi}\). In addition, the Dirac delta function satisfies

$$ \delta(\boldsymbol{x})=\left\{\begin{array}{ll} \infty, & \boldsymbol{x}=\mathbf{0} \\ 0, & \boldsymbol{x} \neq \mathbf{0} \end{array}\right. $$

and when a constant function \(\phi(\boldsymbol{x})=1\) is considered, this function has the following property:

$$ \int_{\Omega} \delta(\boldsymbol{\xi}) \mathrm{d} \boldsymbol{\xi}=1 $$

The \( \delta\) function can be regarded as a very cuspidal function,*1 where the function value becomes \( \infty\) or zero at the origin. In the SPH method, instead of this cuspidal function, a smooth function is used for approximation.

$$ \phi(\boldsymbol{x})=\int_{\Omega} \phi(\boldsymbol{\xi}) \delta(\boldsymbol{\xi}-\boldsymbol{x}) \mathrm{d} \boldsymbol{\xi} $$ $$ \approx \int_{\Omega} \phi(\boldsymbol{\xi}) w_{h}^{\mathrm{SPH}}(|\boldsymbol{\xi}-\boldsymbol{x}|) \mathrm{d} \boldsymbol{\xi}=: \phi^{\mathrm{SPH}}(\boldsymbol{x}) $$

*1 It is a unique function not only in its shape but also in mathematic properties, and is called a generalized function.


\(\phi^{\mathrm{SPH}}(\boldsymbol{x})\) is described as kernel approximation of \(\phi(\boldsymbol{x})\), and this smooth function \(w_{h}^{\mathrm{SPH}}\) is called the “kernel function.” It is assumed that \(h\) indicates a region (radius) to determine the target for smoothing, and in the SPH method, from the functions that satisfy the conditions below, a kernel function is selected.

$$ \text { (Compact support condition) } \quad w_{h}^{\mathrm{SPH}}(|\boldsymbol{\xi}-\boldsymbol{x}|)=0, \quad|\boldsymbol{\xi}-\boldsymbol{x}| \geq h $$ $$ \text { (Unity condition) }\quad \int_{\Omega} w_{h}^{\mathrm{SPH}}(|\boldsymbol{\xi}-\boldsymbol{x}|) d \boldsymbol{\xi}=1 $$

The above approximation Eq. (5) approximates the integral representation Eq. (4) using the delta function, and the compact support condition (6) and the unity condition (7) correspond to property (2) and (3) of the Dirac delta function, respectively. In the original integral representation, only the function value of the observation point \( \boldsymbol{x}\) is extracted from the dummy variable \(\boldsymbol{\xi}\), and therefore, the delta function \(\delta(\boldsymbol{\xi}-\boldsymbol{x}) \) has been utilized. In contrast, in the kernel approximation in the SPH method, points in the neighborhood of the observation point are also incorporated for smoothing, which means that the weighted average is obtained using the values in the neighborhood \(\boldsymbol{\xi}\), where \( |\boldsymbol{\xi}-\boldsymbol{x}|\) denoting the relative distance from the observation point \( \boldsymbol{x}\) is within a fixed radius \(h\). The compact support condition represents that the target for the weighted average is limited to the neighborhood within the fixed radius \(h\). The weight for such smoothing is the kernel function \(w_{h}^{\mathrm{SPH}}(|\boldsymbol{\xi}-\boldsymbol{x}|) \), where weight is given according to relative distance. The unity condition is a necessary condition to converge the kernel function to the delta function when the radius of influence \(h\) is converged to zero. In other words, this means that as the radius of influence is reduced, the weight of each point increases to satisfy the unity condition, and if the radius of influence is converged to the observation point as the limit, it should have infinite weight—in other words, it should be the Dirac delta function.


A kernel approximation in the integral form using this smooth function is the basis of the SPH method, and also became the origin of the term Smoothed Particle Hydrodynamics (SPH). Now, using the Taylor expansion, it is examined whether the equation, in which the Delta function has been changed to a kernel function, sufficiently holds as an approximation. Based on the compact support condition, the dummy variable \(\boldsymbol{\xi}\) appearing in the integrand could be considered as a point in the neighborhood of the observation point \( \boldsymbol{x}\). When Taylor expansion is performed around the observation point, where the function value in the neighborhood of this observation point \(\phi(\boldsymbol{\xi})=\phi(\boldsymbol{x}+\Delta \xi)\) is fixed,*2 the following equations are obtained:

$$ \phi(\boldsymbol{\xi})=\phi(\boldsymbol{x})+\frac{\partial \phi}{\partial \boldsymbol{\xi}}(\boldsymbol{x}) \cdot \Delta \boldsymbol{\xi}+\frac{1}{2} \Delta \boldsymbol{\xi} \cdot \frac{\partial^{2} \phi}{\partial \boldsymbol{\xi}^{2}}(\boldsymbol{x}) \Delta \boldsymbol{\xi}+\mathcal{O}\left(\Delta \boldsymbol{\xi}^{3}\right) $$ $$ =\phi(\boldsymbol{x})+(\boldsymbol{\xi}-\boldsymbol{x}) \cdot \frac{\partial \phi}{\partial \boldsymbol{\xi}}(\boldsymbol{x})+\mathcal{O}\left(\Delta \boldsymbol{\xi}^{2}\right) $$

*2 The same observation point is used at all points in the neighborhood to perform Taylor expansion.

where \(\Delta \boldsymbol{\xi}=\boldsymbol{\xi}-\boldsymbol{x}\) and hereafter, the approximation shall be conducted up to the first-order term of \(\Delta \boldsymbol{\xi}\).

When this Taylor expansion Eq. (9) is substituted for Eq. (5), the following equations are obtained:

$$ \begin{equation} \begin{aligned} \phi^{\mathrm{SPH}}(\boldsymbol{x}) &=\int_{\Omega} \phi(\boldsymbol{\xi}) w_{h}^{\mathrm{SPH}}(|\boldsymbol{\xi}-\boldsymbol{x}|) \mathrm{d} \boldsymbol{\xi} \\ &=\int_{\Omega}\left\{\phi(\boldsymbol{x})+(\boldsymbol{\xi}-\boldsymbol{x}) \cdot \frac{\partial \phi}{\partial \boldsymbol{\xi}}(\boldsymbol{x})\right\} w_{h}^{\mathrm{SPH}}(|\boldsymbol{\xi}-\boldsymbol{x}|) \mathrm{d} \boldsymbol{\xi}+R \end{aligned} \end{equation} $$ $$ =\phi(\boldsymbol{x}) \int_{\Omega} w_{h}^{\mathrm{SPH}}(|\boldsymbol{\xi}-\boldsymbol{x}|) \mathrm{d} \boldsymbol{\xi}+\frac{\partial \phi}{\partial \boldsymbol{\xi}}(\boldsymbol{x}) \cdot \int_{\Omega}(\boldsymbol{\xi}-\boldsymbol{x}) w_{h}^{\mathrm{SPH}}(|\boldsymbol{\xi}-\boldsymbol{x}|) \mathrm{d} \boldsymbol{\xi}+R $$ $$ =\phi(\boldsymbol{x}) \int_{\Omega} w_{h}^{\mathrm{SPH}}\left(\left|\boldsymbol{\xi}^{\prime}\right|\right) \mathrm{d} \boldsymbol{\xi}^{\prime}+\frac{\partial \phi}{\partial \boldsymbol{\xi}}(\boldsymbol{x}) \cdot \int_{\Omega} \boldsymbol{\xi}^{\prime} w_{h}^{\mathrm{SPH}}\left(\left|\boldsymbol{\xi}^{\prime}\right|\right) \mathrm{d} \boldsymbol{\xi}^{\prime}+R $$

\(R\) is an error related to \(\mathcal{O}\left(\Delta \boldsymbol{\xi}^{2}\right)\) in Eq. (9). In the deformation from Eq. (10) to Eq. (11), since the integrand is considered the function of the dummy variable \(\boldsymbol{\xi}\), \(\phi(\boldsymbol{x})\) and \(\frac{\partial \phi}{\partial \boldsymbol{\xi}}(\boldsymbol{x})\) are values at the fixed observation point \(\boldsymbol{x}\), and it is possible to move them outside the integration. Moreover, in Eq. (11) and Eq. (12), \(\boldsymbol{\xi}-\boldsymbol{x}\) denotes a relative position centering on the observation point \(\boldsymbol{x}\) and is expressed by \(\boldsymbol{\xi}^{\prime}\). Ultimately, to ensure the kernel approximation \(\phi^{\mathrm{SPH}}(\boldsymbol{x})\) shown in Eq. (11) holds as an approximation of the function \(\phi(\boldsymbol{x})\), the following two conditions need to be satisfied:

$$ \int_{\Omega} w_{h}^{\mathrm{SPH}}\left(\left|\boldsymbol{\xi}^{\prime}\right|\right) \mathrm{d} \boldsymbol{\xi}^{\prime}=1 $$ $$ \int_{\Omega} \boldsymbol{\xi}^{\prime} w_{h}^{\mathrm{SPH}}\left(\left|\boldsymbol{\xi}^{\prime}\right|\right) \mathrm{d} \boldsymbol{\xi}^{\prime}=0 $$

The first equation indicates that the value of the kernel function defined around the observation point \(\boldsymbol{x}\) should be “1” when it is integrated into the region; this is a condition for the kernel approximation referred to above as the “unity condition.” The second condition is satisfied if the kernel function is an even function (a function in which variables are symmetric with respect to the origin and function values are equal)

$$ w_{h}^{\mathrm{SPH}}\left(\left|\boldsymbol{\xi}^{\prime}\right|\right)=w_{h}^{\mathrm{SPH}}\left(\left|-\boldsymbol{\xi}^{\prime}\right|\right) $$

however, since the kernel function has been set as a function of distance \(r=|\boldsymbol{\xi}-\boldsymbol{x}|\), this property is automatically satisfied. Because of this kernel function property, \(\boldsymbol{\xi}^{\prime} w_{h}^{\mathrm{SPH}}\left(\left|\boldsymbol{\xi}^{\prime}\right|\right)\) becomes an odd function, indicating that the above condition Eq. (14) is satisfied. Lastly, to meet the demand to reduce the computation load as much as possible during a numerical simulation, i.e., by setting a narrow range for region integration, the compact support condition has been added. Now, all necessary conditions (and demands) given to the kernel function have been prepared. All the conditions are summarized below.

(1. Unity condition)
$$ \begin{eqnarray*} \int_{\Omega} w_{h}^{\mathrm{SPH}}(|\boldsymbol{\xi}-\boldsymbol{x}|) \mathrm{d} \boldsymbol{\xi}=1 \end{eqnarray*} $$
(2.Even function property (a function of distance r is assumed))
$$ \begin{eqnarray*} & w_{h}^{\mathrm{SPH}}(|\boldsymbol{\xi}-\boldsymbol{x}|)=w_{h}^{\mathrm{SPH}}(|\boldsymbol{x}-\boldsymbol{\xi}|) \\ & \rightarrow \nabla w_{h}^{\mathrm{SPH}}(|\boldsymbol{\xi}-\boldsymbol{x}|)=-\nabla w_{h}^{\mathrm{SPH}}(|\boldsymbol{x}-\boldsymbol{\xi}|) \end{eqnarray*} $$
(3.Compact support condition)
$$ \begin{eqnarray*} & w_{h}^{\mathrm{SPH}}(|\boldsymbol{\xi}-\boldsymbol{x}|)=0, \quad|\boldsymbol{\xi}-\boldsymbol{x}| \geq h \end{eqnarray*} $$
(4.Non-negativity condition)
$$ \begin{eqnarray*} w_{h}^{\mathrm{SPH}}(|\boldsymbol{\xi}-\boldsymbol{x}|)\verb|>|0, \quad|\boldsymbol{\xi}-\boldsymbol{x}|\verb|<|h \end{eqnarray*} $$
(5.Condition of convergence to the delta function)
$$ \begin{eqnarray*} \lim _{h \rightarrow 0} w_{h}^{\mathrm{SPH}}(|\boldsymbol{\xi}-\boldsymbol{x}|)=\delta(|\boldsymbol{\xi}-\boldsymbol{x}|) \end{eqnarray*} $$
(6.Monotone decreasing condition)
$$ \begin{eqnarray*} \frac{\mathrm{d} w_{h}^{\mathrm{SPH}}}{\mathrm{d} r}\verb|<|0,0\verb|<|r\verb|<|h \end{eqnarray*} $$
(7.Differential continuity)
A function continuous up to the second derivative is often selected.

As the kernel function \(w_{h}^{\mathrm{SPH}}\) that satisfies the above seven conditions, piecewise polynomials called “cubic spline function” and “quintic spline function” are often used in the SPH numerical simulations in accordance with the policy to select the simplest possible functions (see Fig. 1).

$$ w_{h}^{\mathrm{SPHc}}(r)=\alpha_{d}^{\mathrm{c}}\left\{\begin{array}{ll} 1-6(r / h)^{2}+6(r / h)^{3}, & 0 \leq r / h<1 / 2, \\ 2(1-r / h)^{3}, & 1 / 2 \leq r / h<1, \\ 0, & r / h \geq 1, \end{array}\right. $$ $$ w_{h}^{\mathrm{SPHq}}(r)=\alpha_{d}^{\mathrm{q}}\left\{\begin{array}{ll} (1-r / h)^{5}-6(2 / 3-r / h)^{5}+15(1 / 3-r / h)^{5}, & 0 \leq r / h<1 / 3, \\ (1-r / h)^{5}-6(2 / 3-r / h)^{5}, & 1 / 3 \leq r / h<2 / 3, \\ (1-r / h)^{5}, & 2 / 3 \leq r / h<1, \\ 0, & r / h \geq 1 . \end{array}\right. $$
Picture
Fig. 1 Graph of weight functions used in the SPH method (1 is set for h, 2 for d)

Here, \(\alpha_{d}^{c}\) and \(\alpha_{d}^{\mathrm{q}}\) are the constants to satisfy the unity condition (7), and can be determined in advance by setting the radius of influence as follows: \(\alpha_{d}^{c}=10 / 7 \pi h^{2} \) (two dimensions) or \(1 / \pi h^{3} \) (three dimensions); \(\alpha_{d}^{\mathrm{q}}=15309 / 478 \pi h^{2}\) (two dimensions) or \(2187 / 40 \pi h^{3} \) (three dimensions). Since the aforementioned spline functions are piecewise polynomials requiring case division, the Wendland function [8] requiring no case division is also sometimes used.

Reference:
[8] H. Wendland: Piecewise polynomial, positive definite and compactly supported radial functions of minimal degree, Advances in Computational Mathematics volume 4, 389–396, 1995.
Growing the particle method, and its present state
1. Present State of the Particle Method
2. Kernel Approximation and Kernel Function Conditions in the SPH Method (Preparation for Spatial Discretization)
Back to Technical Columns INDEX

 About the author

画像
Kyushu University, Department of Civil Engineering
Mitsuteru Asai
​2003, Ph.D.  (Tohoku University, Department of Civil Engineering.)
2003-2005, Post-Doctoral Researcher (The Ohio State University)
2005-2007, Assistant Professor (Ritsumeikan University)
2007-, Associate Professor (Kyushu University)
(Sitemap)
​Home
Case Examples
 - Particleworks Case Examples
Technical
 - Technical Column
 - References
 - Papers, Lectures
 - Particle Method Glossary
Contact
 - Developer, Partners
 - Contact Us
Privacy Policy
Terms of Use
GDPR PRIVACY POLICY
(Related Sites)
Prometech Software Site
Granuleworks Site
Prometech Simulation
​Conference Site
GDEP Solutions Site
Contact Us
[Developer, Main Domestic / Global Dealer​]
  Prometech Software, Inc.
Prometech Software, Inc.
​URL: www.prometech.co.jp
​
E-mail: web@prometech.co.jp
(Sitemap)
Home
事例
 - 解析事例
Learning
 - 粒子法・MPS法
​ - 技術コラム
 - 粒子法用語集
 - 参考文献・ウェブサイト
 - 論文・講演
お問い合わせ・ご相談
 - 導入の流れとライセンス形態
 - Particleworks / Granuleworks
   プリインストールGPU搭載
   ワークステーション
 - 開発元・パートナー
 - お問い合わせ
プライバシーポリシー
利用規約
GDPR プライバシーポリシー
(Related Sites)
Prometech Software サイト
​Granuleworksサイト
​Prometech Simulation
Conference サイト
Particleworks Europe サイト
プロメテックCGリサーチ サイト
GDEP Solutions サイト
-  動作確認済み GPU搭載ワークステーション
HPC WORLD サイト
お問い合わせフォーム
[開発元・国内、海外総販売店]
  プロメテック・ソフトウェア株式会社
Prometech Software, Inc.
URL: www.prometech.co.jp
E-mail: web@prometech.co.jp

ⓒPrometech Software, Inc.
  • Home
  • Home
  • Case Examples
    • Particleworks解析事例
    • Particleworks Case Examples
  • Technical
    • 粒子法・MPS法
    • 技術コラム >
      • DX時代の製品開発プロセスとCAEの重要性 >
        • 第1回 序 略歴とコラム紹介
        • 第2回 DXとデジタルエンジニアリング
        • 第3回 製品開発プロセスの目指す姿
        • 第4回 DX時代のCAE
        • 第5回 評価CAEの概要と課題
        • 第6回 評価CAEの課題解決手法
        • 第7回 企画CAEの概要と課題
        • 第8回 企画CAEの運用と応用
        • 第9回 設計CAEの概要と課題
        • 第10回 設計CAEの課題解決の進め方
        • 第11回 開発プロセス運用の仕組み作り
        • 第12回 まとめと変革の時代に求められるエンジニア像
      • 粒子法のいま、そして未来へ >
        • 第1回 粒子法のいま
        • 第2回 SPH法におけるカーネル近似とカーネル関数の条件
        • 第3回 SPH法における空間離散化
      • 粒子法の非圧縮条件とは
      • 粒子法入門 >
        • 第1回 粒子法って何?
        • 第2回 粒子法は、他の方法とどう違うか
        • 第3回 粒子法の大きさと質量について
        • ​第4回 「粒子の動かし方」と「加速度の求め方」について
        • ​第5回 計算時間を短縮する方法について
    • Technical Column >
      • Growing the particle method, and its present state >
        • 1. Present State of the Particle Method
        • 2. Kernel Approximation and Kernel Function Conditions in the SPH Method (Preparation for Spatial Discretization)
      • Incompressibility of the particle method
      • Introduction to the particle method >
        • 1. What is a particle method?
        • 2. In what ways is the particle method different from other methods?
        • 3. Mass and volume of particles
        • 4. How to move particles and how to calculate accelerations of particles
        • 5. How to shorten the simulation time
    • 粒子法用語集
    • Particle Method Glossary
    • 参考文献・ウェブサイト
    • Reference Book/URL
    • 論文・講演
  • Contact
    • 導入の流れとライセンス形態
    • Particleworks / GranuleworksプリインストールGPU搭載ワークステーション
    • 開発元・パートナー
    • Developers, Partners
    • お問い合わせ
    • Contact Us