凸多边形/多面体上线性函数平均值的分解(积分点和积分权重)
背景
在一般凸多边形或者多面体单元的二阶精度的有限体积格式,例如MUSCL,或者DG格式,需要在单元边界的中点上取数值通量,权重用于稳定性和保正性的证明。下面是给出一般凸多胞体(多边形或者多面体)的积分点和权重的推导。
结论
设在单元(控制体)
V
i
V_i
Vi 上守恒量
U
U
U的平均值为
U
i
U_i
Ui ,在单元
V
i
V_i
Vi 上重构的线性函数为,
U
(
x
)
=
U
i
+
∇
U
⋅
(
x
−
x
i
)
U(x) = U_i + \nabla U \cdot(\boldsymbol{x}-\boldsymbol{x}_i)
U(x)=Ui+∇U⋅(x−xi)
其中
x
i
\boldsymbol{x}_i
xi 是控制体
V
i
V_i
Vi 的重心。控制点
x
i
j
x_{ij}
xij 为控制体
V
i
V_i
Vi 的边界
S
i
j
S_{ij}
Sij 的重心,下标
j
∈
N
(
i
)
j \in N(i)
j∈N(i) ,
N
(
i
)
N(i)
N(i) 表示单元
V
i
V_i
Vi 的邻居指标集合。控制点的权重
w
i
j
w_{ij}
wij 满足,
∑
j
∈
N
(
i
)
w
i
j
=
1
(1)
\sum_{j \in N(i)} w_{ij}=1 \tag{1}
j∈N(i)∑wij=1(1)
其中权重
w
i
j
=
(
x
i
j
−
x
i
)
⋅
n
i
j
∣
S
i
j
∣
D
∣
V
i
∣
(2)
w_{ij} = \frac{(\boldsymbol{x}_{ij}-\boldsymbol{x}_i)\cdot\boldsymbol{n}_{ij}|S_{ij}|}{D|V_i|} \tag{2}
wij=D∣Vi∣(xij−xi)⋅nij∣Sij∣(2)
其中,
n
i
j
\boldsymbol{n}_{ij}
nij 是边界
S
i
j
S_{ij}
Sij 指向单元外部的单位法向量,
∣
S
i
j
∣
|S_{ij}|
∣Sij∣ 是边界
S
i
j
S_{ij}
Sij 的长度或者面积。在几何上
(
x
i
j
−
x
i
)
⋅
n
i
j
(\boldsymbol{x}_{ij}-\boldsymbol{x}_i)\cdot\boldsymbol{n}_{ij}
(xij−xi)⋅nij 表示重心
x
i
x_i
xi 到
S
i
j
S_{ij}
Sij 所在直线或者平面上的距离。
线性函数
U
(
x
)
U(x)
U(x) 的数值积分满足,
U
i
=
∑
j
w
i
j
U
i
j
(3)
U_i = \sum_{j} w_{ij}U_{ij} \tag{3}
Ui=j∑wijUij(3)
其中,
U
i
j
=
U
(
x
i
j
)
U_{ij}=U(\boldsymbol{x}_{ij})
Uij=U(xij) 。
证明
预备知识
Gauss-Green 定理
对于向量函数
F
∈
R
n
\boldsymbol{F}\in \mathbb{R}^n
F∈Rn ,在区域
Ω
\Omega
Ω 上,有Gauss-Green公式,
∫
Ω
∇
⋅
F
d
x
=
∫
∂
Ω
F
⋅
n
d
S
(
x
)
(4)
\int_{\Omega} \nabla \cdot \boldsymbol{F} d\boldsymbol{x} =\int_{\partial \Omega} \boldsymbol{F} \cdot \boldsymbol{n} dS(\boldsymbol{x}) \tag{4}
∫Ω∇⋅Fdx=∫∂ΩF⋅ndS(x)(4)
Gauss-Green 变形公式
对于标量函数
u
,
v
u,v
u,v ,变形公式为
∫
Ω
∇
u
v
d
x
=
−
∫
Ω
u
∇
v
d
x
+
∫
∂
Ω
u
v
n
d
x
(5)
\int_{\Omega} \nabla u v d \boldsymbol{x} = -\int_{\Omega} u \nabla v d \boldsymbol{x} + \int_{\partial \Omega} u v \boldsymbol{n} d \boldsymbol{x} \tag{5}
∫Ω∇uvdx=−∫Ωu∇vdx+∫∂Ωuvndx(5)
证明结论(1)(2)
多胞体的重心公式满足
x
i
=
1
∣
V
i
∣
∫
V
i
x
d
V
(6)
\boldsymbol{x}_i = \frac{1}{|V_i|} \int_{V_{i}} \boldsymbol{x} d V \tag{6}
xi=∣Vi∣1∫VixdV(6)
控制体
V
i
V_i
Vi的体积公式为
∣
V
i
∣
=
∫
V
i
1
d
V
=
∫
V
i
∇
⋅
x
D
d
V
=
1
D
∫
∂
V
i
x
⋅
n
d
S
=
1
D
∑
j
x
i
j
⋅
n
i
j
∣
S
i
j
∣
(7)
|V_i|= \int_{V_i} 1 dV = \int_{V_i} \frac{\nabla\cdot \boldsymbol{x}}{D}dV = \frac{1}{D}\int_{\partial V_i} \boldsymbol{x} \cdot \boldsymbol{n} dS = \frac{1}{D} \sum_j \boldsymbol{x}_{ij} \cdot \boldsymbol{n}_{ij}|S_{ij}| \tag{7}
∣Vi∣=∫Vi1dV=∫ViD∇⋅xdV=D1∫∂Vix⋅ndS=D1j∑xij⋅nij∣Sij∣(7)
这里倒数第二个等号使用了Gauss-Green公式(4),导数第一个等号成立的是因为边界
S
i
j
S_{ij}
Sij是线性函数,在边界上取重心点就可以求出在整个边界上的积分,后面在单元边界
∂
V
i
\partial V_i
∂Vi上离散时都会使用到,下面再赘述。
单元的重心
x
i
\boldsymbol{x}_i
xi 是一个常数,因而
∇
⋅
x
i
=
0
\nabla \cdot \boldsymbol{x}_i =0
∇⋅xi=0 , 再套用公式(7)可以得到
0
=
∫
V
i
0
d
V
=
1
D
∫
V
i
∇
⋅
x
i
d
V
=
1
D
∑
j
x
i
⋅
n
i
j
∣
S
i
j
∣
(8)
0 = \int_{V_i}0 d V = \frac{1}{D}\int_{V_i} \nabla\cdot\boldsymbol{x}_i d V = \frac{1}{D} \sum_j \boldsymbol{x}_i \cdot \boldsymbol{n}_{ij}|S_{ij}| \tag{8}
0=∫Vi0dV=D1∫Vi∇⋅xidV=D1j∑xi⋅nij∣Sij∣(8)
将公式(7)(8)作差得到
∣
V
i
∣
=
1
D
∑
j
x
i
j
⋅
n
i
j
∣
S
i
j
∣
−
1
D
∑
j
x
i
⋅
n
i
j
∣
S
i
j
∣
=
1
D
∑
j
(
x
i
j
−
x
i
)
⋅
n
i
j
∣
S
i
j
∣
|V_i| = \frac{1}{D} \sum_j \boldsymbol{x}_{ij} \cdot \boldsymbol{n}_{ij}|S_{ij}| - \frac{1}{D} \sum_j \boldsymbol{x}_i \cdot \boldsymbol{n}_{ij}|S_{ij}| = \frac{1}{D} \sum_j (\boldsymbol{x}_{ij}-\boldsymbol{x}_i) \cdot \boldsymbol{n}_{ij}|S_{ij}|
∣Vi∣=D1j∑xij⋅nij∣Sij∣−D1j∑xi⋅nij∣Sij∣=D1j∑(xij−xi)⋅nij∣Sij∣
由此可以得到结论(1)(2).
证明结论(3)
对于标量函数
U
U
U满足
U
=
1
D
D
U
=
1
D
(
∇
⋅
x
)
U
=
1
D
[
∇
⋅
(
U
x
)
−
x
⋅
∇
U
]
(9)
U= \frac{1}{D} D U = \frac{1}{D}(\nabla \cdot \boldsymbol{x})U = \frac{1}{D}[\nabla\cdot(U\boldsymbol{x})-\boldsymbol{x}\cdot \nabla U] \tag{9}
U=D1DU=D1(∇⋅x)U=D1[∇⋅(Ux)−x⋅∇U](9)
单元的平均值
U
i
U_i
Ui 的定义为
U
i
=
1
∣
V
i
∣
∫
V
i
U
d
V
U_i = \frac{1}{|V_i|}\int_{V_i} U d V
Ui=∣Vi∣1∫ViUdV
代入公式(9)得到
U
i
=
1
D
∣
V
i
∣
∫
V
i
[
∇
⋅
(
U
x
)
−
x
⋅
∇
U
]
d
V
=
1
D
∣
V
i
∣
[
∫
V
i
∇
⋅
(
U
x
)
d
V
−
∇
U
⋅
∫
V
i
x
d
V
]
(10)
U_i =\frac{1}{D|V_i|}\int_{V_i} [\nabla\cdot(U\boldsymbol{x})-\boldsymbol{x}\cdot \nabla U] dV =\frac{1}{D|V_i|}[\int_{V_i} \nabla\cdot(U\boldsymbol{x})dV-\nabla U \cdot \int_{V_i} \boldsymbol{x} dV ] \tag{10}
Ui=D∣Vi∣1∫Vi[∇⋅(Ux)−x⋅∇U]dV=D∣Vi∣1[∫Vi∇⋅(Ux)dV−∇U⋅∫VixdV](10)
这里第二项
∇
U
\nabla U
∇U 可以放在积分号之外是因为
∇
U
\nabla U
∇U 是一个常数。因此对于梯度
∇
U
\nabla U
∇U同样有
∇
U
=
1
∣
V
i
∣
∫
V
i
∇
U
d
V
=
1
V
i
∫
∂
V
i
U
n
d
A
=
1
∣
V
i
∣
∑
j
U
i
j
n
i
j
∣
S
i
j
∣
(11)
\nabla U =\frac{1}{|V_i|}\int_{V_i} \nabla U dV = \frac{1}{V_i} \int_{\partial V_i} U \boldsymbol{n} dA =\frac{1}{|V_i|} \sum_j U_{ij}\boldsymbol{n}_{ij}|S_{ij}| \tag{11}
∇U=∣Vi∣1∫Vi∇UdV=Vi1∫∂ViUndA=∣Vi∣1j∑Uijnij∣Sij∣(11)
这里第二个等号使用到了Gauss-Green公式的变形公式(5).
公式(10)的第一项应用Gauss-Green公式(4),第二项代入重心公式(6),可以得到
U
i
=
1
D
∣
V
i
∣
[
∫
∂
V
i
U
x
⋅
n
d
S
−
∇
U
⋅
x
i
]
U_i = \frac{1}{D|V_i|}\left [ \int_{\partial V_i} U \boldsymbol{x} \cdot \boldsymbol{n} dS -\nabla U \cdot \boldsymbol{x}_i\right]
Ui=D∣Vi∣1[∫∂ViUx⋅ndS−∇U⋅xi]
代入公式(11)可以得到
U
i
=
1
D
∣
V
i
∣
[
∫
∂
V
i
U
x
⋅
n
d
S
−
x
i
⋅
∫
V
i
∇
U
d
V
]
U_i = \frac{1}{D|V_i|}\left [ \int_{\partial V_i} U \boldsymbol{x} \cdot \boldsymbol{n} dS -\boldsymbol{x}_i \cdot \int_{V_i} \nabla U dV\right]
Ui=D∣Vi∣1[∫∂ViUx⋅ndS−xi⋅∫Vi∇UdV]
离散之后得到
U
i
=
1
D
∣
V
i
∣
[
∑
j
U
i
j
(
x
i
j
⋅
n
i
j
)
∣
S
i
j
∣
−
x
i
⋅
∑
j
U
i
j
∣
S
i
j
∣
n
i
j
]
=
1
D
∣
V
i
∣
∑
j
(
x
i
j
−
x
i
)
⋅
n
i
j
∣
S
i
j
∣
U
i
j
=
∑
j
w
i
j
U
i
j
U_i = \frac{1}{D|V_i|} \left[ \sum_j U_{ij}(\boldsymbol{x}_{ij}\cdot \boldsymbol{n}_{ij})|S_{ij}|-\boldsymbol{x}_i \cdot \sum_j U_{ij} |S_{ij}|\boldsymbol{n}_{ij} \right] = \frac{1}{D|V_i|} \sum_j (\boldsymbol{x}_{ij}-\boldsymbol{x}_{i})\cdot \boldsymbol{n}_{ij} |S_{ij}|U_{ij}=\sum_j w_{ij}U_{ij}
Ui=D∣Vi∣1[j∑Uij(xij⋅nij)∣Sij∣−xi⋅j∑Uij∣Sij∣nij]=D∣Vi∣1j∑(xij−xi)⋅nij∣Sij∣Uij=j∑wijUij
得证结论
(
3
)
(3)
(3).