聚合超椭球（Poly-super-ellipsoid） [This page is prepared in Chinese, and will be part of a book aiming at beginners of DEM.]
聚合超椭球（Poly-super-ellipsoid）
超椭球的表面可表示为 (Zhao & Zhao, 2019)
f ( x , y , z ) = ( ∣ x r x ∣ 2 ϵ 1 + ∣ y r y ∣ 2 ϵ 1 ) ϵ 1 ϵ 2 + ∣ z r z ∣ 2 ϵ 2 − 1 = 0 (1) \begin{aligned}
f(x, y, z) = \left( \left| \frac{x}{r_{x}} \right|^{\frac{2}{\epsilon_{1}}}
+ \left| \frac{y}{r_{y}} \right|^{\frac{2}{\epsilon_{1}}} \right)^{\frac{\epsilon_{1}}{\epsilon_{2}}}
+ \left| \frac{z}{r_{z}} \right|^{\frac{2}{\epsilon_{2}}} - 1 = 0
\tag{1}
\end{aligned} f ( x , y , z ) = ( r x x ϵ 1 2 + r y y ϵ 1 2 ) ϵ 2 ϵ 1 + r z z ϵ 2 2 − 1 = 0 ( 1 ) 其中，r x r_{x} r x r y r_{y} r y r z r_{z} r z x x x y y y z z z ϵ 1 \epsilon_{1} ϵ 1 ϵ 2 \epsilon_{2} ϵ 2 ϵ 1 \epsilon_{1} ϵ 1 ϵ 2 \epsilon_{2} ϵ 2 r x r_{x} r x r x − r_{x^-} r x − r x + r_{x^+} r x + r y r_{y} r y r z r_{z} r z
支撑映射（Support mapping）
参考笛卡尔坐标系和球坐标系之间的坐标变换，超椭球的球坐标参数化可表示为（注，以第一象限为例,正负号已省略）
x = r x ( sin θ ) ϵ 2 ( cos φ ) ϵ 1 y = r y ( sin θ ) ϵ 2 ( sin φ ) ϵ 1 z = r z ( cos θ ) ϵ 2 \begin{aligned}
x &= r_x (\sin\theta)^{\epsilon_2} (\cos\varphi)^{\epsilon_1} \\
y &= r_y (\sin\theta)^{\epsilon_2} (\sin\varphi)^{\epsilon_1} \\
z &= r_z (\cos\theta)^{\epsilon_2}
\end{aligned} x y z = r x ( sin θ ) ϵ 2 ( cos φ ) ϵ 1 = r y ( sin θ ) ϵ 2 ( sin φ ) ϵ 1 = r z ( cos θ ) ϵ 2 其中，θ \theta θ φ \varphi φ x x x y y y z z z θ \theta θ φ \varphi φ
x θ = ϵ 2 r x ( sin θ ) ϵ 2 − 1 ( cos θ ) ( cos φ ) ϵ 1 y θ = ϵ 2 r y ( sin θ ) ϵ 2 − 1 ( cos θ ) ( sin φ ) ϵ 1 z θ = ϵ 2 r z ( cos θ ) ϵ 2 − 1 ( − sin θ ) \begin{aligned}
x_\theta &= \epsilon_2 r_x (\sin\theta)^{\epsilon_2-1} (\cos\theta)
(\cos\varphi)^{\epsilon_1} \\
y_\theta &= \epsilon_2 r_y (\sin\theta)^{\epsilon_2-1} (\cos\theta)
(\sin\varphi)^{\epsilon_1} \\
z_\theta &= \epsilon_2 r_z (\cos\theta)^{\epsilon_2-1} (-\sin\theta)
\end{aligned} x θ y θ z θ = ϵ 2 r x ( sin θ ) ϵ 2 − 1 ( cos θ ) ( cos φ ) ϵ 1 = ϵ 2 r y ( sin θ ) ϵ 2 − 1 ( cos θ ) ( sin φ ) ϵ 1 = ϵ 2 r z ( cos θ ) ϵ 2 − 1 ( − sin θ ) 以及
x φ = ϵ 1 r x ( sin θ ) ϵ 2 ( cos φ ) ϵ 1 − 1 ( − sin φ ) y φ = ϵ 1 r y ( sin θ ) ϵ 2 ( sin φ ) ϵ 1 − 1 ( cos φ ) z φ = 0 \begin{aligned}
x_\varphi &= \epsilon_1 r_x (\sin\theta)^{\epsilon_2}
(\cos\varphi)^{\epsilon_1-1} (-\sin\varphi) \\
y_\varphi &= \epsilon_1 r_y (\sin\theta)^{\epsilon_2}
(\sin\varphi)^{\epsilon_1-1} (\cos\varphi) \\
z_\varphi &= 0
\end{aligned} x φ y φ z φ = ϵ 1 r x ( sin θ ) ϵ 2 ( cos φ ) ϵ 1 − 1 ( − sin φ ) = ϵ 1 r y ( sin θ ) ϵ 2 ( sin φ ) ϵ 1 − 1 ( cos φ ) = 0 求两个方向向量的叉乘,并化简，可得到超椭球表面的法方向为（注，法方向各分量同时除于了 ϵ 1 ϵ 2 r x r y r z ( sin θ ) 2 ϵ 2 − 1 ( cos θ ) ϵ 2 − 1 ( sin φ ) ϵ 1 − 1 ( cos φ ) ϵ 1 − 1 \epsilon_1 \epsilon_2 r_x r_y r_z (\sin\theta)^{2\epsilon_2-1} (\cos\theta)^{\epsilon_2-1} (\sin\varphi)^{\epsilon_1-1} (\cos\varphi)^{\epsilon_1-1} ϵ 1 ϵ 2 r x r y r z ( sin θ ) 2 ϵ 2 − 1 ( cos θ ) ϵ 2 − 1 ( sin φ ) ϵ 1 − 1 ( cos φ ) ϵ 1 − 1
n x = k r x ( sin θ ) 2 − ϵ 2 ( cos φ ) 2 − ϵ 1 n y = k r y ( sin θ ) 2 − ϵ 2 ( sin φ ) 2 − ϵ 1 n z = k r z ( cos θ ) 2 − ϵ 2 \begin{aligned}
n_x &= \frac{k}{r_x} (\sin\theta)^{2-\epsilon_2} (\cos\varphi)^{2-\epsilon_1} \\
n_y &= \frac{k}{r_y} (\sin\theta)^{2-\epsilon_2} (\sin\varphi)^{2-\epsilon_1} \\
n_z &= \frac{k}{r_z} (\cos\theta)^{2-\epsilon_2}
\end{aligned} n x n y n z = r x k ( sin θ ) 2 − ϵ 2 ( cos φ ) 2 − ϵ 1 = r y k ( sin θ ) 2 − ϵ 2 ( sin φ ) 2 − ϵ 1 = r z k ( cos θ ) 2 − ϵ 2 其中，n x n_x n x n y n_y n y n z n_z n z k k k n x n_x n x n y n_y n y n z n_z n z θ \theta θ φ \varphi φ
φ = atan2 ( ( r y n y ) 1 2 − ϵ 1 , ( r x n x ) 1 2 − ϵ 1 ) θ = atan2 ( ( r x n x / ( cos φ ) 2 − ϵ 1 ) 1 2 − ϵ 2 , ( r z n z ) 1 2 − ϵ 2 ) or θ = atan2 ( ( r y n y / ( sin φ ) 2 − ϵ 1 ) 1 2 − ϵ 2 , ( r z n z ) 1 2 − ϵ 2 ) \begin{aligned}
\varphi &= \text{atan2}
\left(
(r_y n_y)^{\frac{1}{2-\epsilon_1}}, (r_x n_x)^{\frac{1}{2-\epsilon_1}}
\right) \\
\theta &= \text{atan2}
\left(
\left(r_x n_x / (\cos\varphi)^{2-\epsilon_1} \right)^{\frac{1}{2-\epsilon_2}},
(r_z n_z)^{\frac{1}{2-\epsilon_2}}
\right) \\
\text{or} \quad
\theta &= \text{atan2}
\left(
\left(r_y n_y / (\sin\varphi)^{2-\epsilon_1} \right)^{\frac{1}{2-\epsilon_2}},
(r_z n_z)^{\frac{1}{2-\epsilon_2}}
\right)
\end{aligned} φ θ or θ = atan2 ( ( r y n y ) 2 − ϵ 1 1 , ( r x n x ) 2 − ϵ 1 1 ) = atan2 ( ( r x n x / ( cos φ ) 2 − ϵ 1 ) 2 − ϵ 2 1 , ( r z n z ) 2 − ϵ 2 1 ) = atan2 ( ( r y n y / ( sin φ ) 2 − ϵ 1 ) 2 − ϵ 2 1 , ( r z n z ) 2 − ϵ 2 1 ) 有向距离场（Signed distance field）
给定查询节点 P P P
x ⃗ Q = ( x Q , y Q , z Q ) = ( c x P , c y P , c z P ) \begin{aligned}
\vec{x}_Q = (x_Q,y_Q,z_Q) = (c x_P,c y_P,c z_P)
\end{aligned} x Q = ( x Q , y Q , z Q ) = ( c x P , c y P , c z P ) 其中，x x x y y y z z z P P P Q Q Q P P P Q Q Q c c c Q Q Q ( 1 ) (1) ( 1 ) Q Q Q ( 1 ) (1) ( 1 ) c c c
c = ( 1 f ( x ⃗ P ) + 1 ) ϵ 2 2 \begin{aligned}
c = \left( \frac{1}{f(\vec{x}_P)+1} \right)^{\frac{\epsilon_2}{2}}
\end{aligned} c = ( f ( x P ) + 1 1 ) 2 ϵ 2 有向距离场可定义为点 P P P Q Q Q
Φ ( x ⃗ P ) = ( c − 1 ) ∣ ∣ x ⃗ P ∣ ∣ \begin{aligned}
\Phi(\vec{x}_P) = (c-1)||\vec{x}_P||
\end{aligned} Φ ( x P ) = ( c − 1 ) ∣∣ x P ∣∣ 其中，∣ ∣ x ⃗ P ∣ ∣ ||\vec{x}_P|| ∣∣ x P ∣∣ x ⃗ P \vec{x}_P x P
聚合超二次曲面（Poly-super-quadrics)
超二次曲面可表示为
f ( x , y , z ) = ∣ x r x ∣ 2 ϵ x + ∣ y r y ∣ 2 ϵ y + ∣ z r z ∣ 2 ϵ z − 1 = 0 (2) \begin{aligned}
f(x, y, z) = \left| \frac{x}{r_{x}} \right|^{\frac{2}{\epsilon_x}}
+ \left| \frac{y}{r_{y}} \right|^{\frac{2}{\epsilon_y}}
+ \left| \frac{z}{r_{z}} \right|^{\frac{2}{\epsilon_z}} - 1 = 0
\tag{2}
\end{aligned} f ( x , y , z ) = r x x ϵ x 2 + r y y ϵ y 2 + r z z ϵ z 2 − 1 = 0 ( 2 ) 其中，r x r_{x} r x r y r_{y} r y r z r_{z} r z x x x y y y z z z ϵ x \epsilon_{x} ϵ x ϵ y \epsilon_{y} ϵ y ϵ z \epsilon_{z} ϵ z ϵ x \epsilon_{x} ϵ x ϵ y \epsilon_{y} ϵ y ϵ z \epsilon_{z} ϵ z r x r_{x} r x r x − r_{x^-} r x − r x + r_{x^+} r x + ϵ x \epsilon_{x} ϵ x ϵ x − \epsilon_{x^-} ϵ x − ϵ x + \epsilon_{x^+} ϵ x + r y r_{y} r y r z r_{z} r z ϵ y \epsilon_{y} ϵ y ϵ z \epsilon_{z} ϵ z
支撑映射（Support mapping）
参考笛卡尔坐标系和球坐标系之间的坐标变换，超椭球的球坐标参数化可表示为
x = r x ( sin θ ) ϵ x ( cos φ ) ϵ x y = r y ( sin θ ) ϵ y ( sin φ ) ϵ y z = r z ( cos θ ) ϵ z \begin{aligned}
x &= r_x (\sin\theta)^{\epsilon_x} (\cos\varphi)^{\epsilon_x} \\
y &= r_y (\sin\theta)^{\epsilon_y} (\sin\varphi)^{\epsilon_y} \\
z &= r_z (\cos\theta)^{\epsilon_z}
\end{aligned} x y z = r x ( sin θ ) ϵ x ( cos φ ) ϵ x = r y ( sin θ ) ϵ y ( sin φ ) ϵ y = r z ( cos θ ) ϵ z 其中，θ \theta θ φ \varphi φ x x x y y y z z z θ \theta θ φ \varphi φ
x θ = ϵ x r x ( sin θ ) ϵ x − 1 ( cos θ ) ( cos φ ) ϵ x y θ = ϵ y r y ( sin θ ) ϵ y − 1 ( cos θ ) ( sin φ ) ϵ y z θ = ϵ z r z ( cos θ ) ϵ z − 1 ( − sin θ ) \begin{aligned}
x_\theta &= \epsilon_x r_x (\sin\theta)^{\epsilon_x-1} (\cos\theta)
(\cos\varphi)^{\epsilon_x} \\
y_\theta &= \epsilon_y r_y (\sin\theta)^{\epsilon_y-1} (\cos\theta)
(\sin\varphi)^{\epsilon_y} \\
z_\theta &= \epsilon_z r_z (\cos\theta)^{\epsilon_z-1} (-\sin\theta)
\end{aligned} x θ y θ z θ = ϵ x r x ( sin θ ) ϵ x − 1 ( cos θ ) ( cos φ ) ϵ x = ϵ y r y ( sin θ ) ϵ y − 1 ( cos θ ) ( sin φ ) ϵ y = ϵ z r z ( cos θ ) ϵ z − 1 ( − sin θ ) 以及
x φ = ϵ x r x ( sin θ ) ϵ x ( cos φ ) ϵ x − 1 ( − sin φ ) y φ = ϵ y r y ( sin θ ) ϵ y ( sin φ ) ϵ y − 1 ( cos φ ) z φ = 0 \begin{aligned}
x_\varphi &= \epsilon_x r_x (\sin\theta)^{\epsilon_x}
(\cos\varphi)^{\epsilon_x-1} (-\sin\varphi) \\
y_\varphi &= \epsilon_y r_y (\sin\theta)^{\epsilon_y}
(\sin\varphi)^{\epsilon_y-1} (\cos\varphi) \\
z_\varphi &= 0
\end{aligned} x φ y φ z φ = ϵ x r x ( sin θ ) ϵ x ( cos φ ) ϵ x − 1 ( − sin φ ) = ϵ y r y ( sin θ ) ϵ y ( sin φ ) ϵ y − 1 ( cos φ ) = 0 求两个方向向量的叉乘,并化简，可得到超椭球表面的法方向为（注，法方向各分量同时除于了 ϵ x ϵ y ϵ z r x r y r z ( sin θ ) ϵ x + ϵ y − 1 ( cos θ ) ϵ z − 1 ( sin φ ) ϵ y − 1 ( cos φ ) ϵ x − 1 \epsilon_x \epsilon_y \epsilon_z r_x r_y r_z (\sin\theta)^{\epsilon_x+\epsilon_y-1} (\cos\theta)^{\epsilon_z-1} (\sin\varphi)^{\epsilon_y-1} (\cos\varphi)^{\epsilon_x-1} ϵ x ϵ y ϵ z r x r y r z ( sin θ ) ϵ x + ϵ y − 1 ( cos θ ) ϵ z − 1 ( sin φ ) ϵ y − 1 ( cos φ ) ϵ x − 1
n x = k ϵ x r x ( sin θ ) 2 − ϵ x ( cos φ ) 2 − ϵ x n y = k ϵ y r y ( sin θ ) 2 − ϵ y ( sin φ ) 2 − ϵ y n z = k ϵ z r z ( cos θ ) 2 − ϵ z \begin{aligned}
n_x &= \frac{k}{\epsilon_x r_x}
(\sin\theta)^{2-\epsilon_x} (\cos\varphi)^{2-\epsilon_x} \\
n_y &= \frac{k}{\epsilon_y r_y}
(\sin\theta)^{2-\epsilon_y} (\sin\varphi)^{2-\epsilon_y} \\
n_z &= \frac{k}{\epsilon_z r_z}
(\cos\theta)^{2-\epsilon_z}
\end{aligned} n x n y n z = ϵ x r x k ( sin θ ) 2 − ϵ x ( cos φ ) 2 − ϵ x = ϵ y r y k ( sin θ ) 2 − ϵ y ( sin φ ) 2 − ϵ y = ϵ z r z k ( cos θ ) 2 − ϵ z 其中，n x n_x n x n y n_y n y n z n_z n z k k k n x n_x n x n y n_y n y n z n_z n z
x = r x ( ϵ x r x n x k ) ϵ x 2 − ϵ x y = r y ( ϵ y r y n y k ) ϵ y 2 − ϵ y z = r z ( ϵ z r z n z k ) ϵ z 2 − ϵ z \begin{aligned}
x &= r_x \left( \frac{\epsilon_x r_x n_x}{k} \right)^{\frac{\epsilon_x}{2-\epsilon_x}} \\
y &= r_y \left( \frac{\epsilon_y r_y n_y}{k} \right)^{\frac{\epsilon_y}{2-\epsilon_y}} \\
z &= r_z \left( \frac{\epsilon_z r_z n_z}{k} \right)^{\frac{\epsilon_z}{2-\epsilon_z}} \\
\end{aligned} x y z = r x ( k ϵ x r x n x ) 2 − ϵ x ϵ x = r y ( k ϵ y r y n y ) 2 − ϵ y ϵ y = r z ( k ϵ z r z n z ) 2 − ϵ z ϵ z 将上式带入公式 ( 2 ) (2) ( 2 ) k k k x x x y y y z z z
Zhao, S., & Zhao, J. (2019). A poly-superellipsoid-based approach on particle morphology for DEM modeling of granular media. International Journal for Numerical and Analytical Methods in Geomechanics , 43 (13), 2147–2169.