Superconvergent and divergence-free finite element methods for the Stokes equation are developed. The velocity and pressure are discretized using $H(\mathrm{div})$-conforming vector elements and discontinuous piecewise polynomials. The discrete formulation employs a weak deviatoric gradient operator built with tangential-normal continuous finite elements for traceless tensors, requiring no stabilization. Optimal and superconvergent error estimates are established. The method connects to nonconforming virtual element and pseudostress-velocity-pressure mixed formulations. Numerical experiments verify the theory.
논문 ID : 2510.14192제목 : Superconvergent and Divergence-Free Finite Element Methods for Stokes Equation저자 : Long Chen, Xuehai Huang, Chao Zhang, Xinyue Zhao분류 : math.NA cs.NA제출 시간 : 2025년 10월 16일 arXiv 제출논문 링크 : https://arxiv.org/abs/2510.14192 본 논문은 Stokes 방정식을 풀기 위한 초수렴 및 무발산 유한요소법을 개발했습니다. 속도장은 H(div)-적합 벡터 요소로 이산화되고, 압력장은 불연속 분할 다항식으로 이산화됩니다. 이산 형식은 무적 텐서 기반의 접선-법선 연속 유한요소로 구성된 약한 편차 기울기 연산자를 사용하며, 안정화 처리가 필요하지 않습니다. 최적 및 초수렴 오차 추정이 확립되었습니다. 본 방법은 비적합 가상요소법 및 의사응력-속도-압력 혼합 형식과 관련이 있습니다. 수치 실험이 이론적 결과를 검증합니다.
Stokes 방정식은 비압축성 유체 운동을 설명하는 기본 방정식으로, 전산유체역학에서 핵심적 위치를 차지합니다. 이 방정식의 수치 해석은 두 가지 주요 과제에 직면합니다:
무발산 제약 : 속도장은 div u = 0을 만족해야 하며, 이는 질량 보존을 나타냅니다압력 견고성 : 수치 방법은 압력항 처리에 견고성을 가져야 합니다고전적인 안정화 유한요소 쌍(예: Taylor-Hood 요소, MINI 요소, 비적합 P1-P0 요소)은 다음과 같은 문제가 있습니다:
점별 무발산 성질 부재 : 무발산 제약을 엄격히 만족할 수 없어 질량 보존성이 약화됩니다압력 견고성 부족 : 오차 추정에서 압력 견고성이 부족합니다수렴 차수 제한 : 일부 방법에서 수렴 차수 손실이 발생합니다본 논문은 엄격한 무발산 성질을 유지하면서 초수렴을 달성할 수 있는 유한요소법을 구성하는 것을 목표로 하며, 동시에 기존 방법에서 필요한 안정화 처리를 피합니다.
새로운 혼합 유한요소 형식 제안 : H(div)-적합 요소로 속도를 이산화하고, 무연속 분할 다항식으로 압력을 이산화합니다약한 편차 기울기 연산자 구성 : 무적 텐서 기반의 접선-법선 연속 유한요소를 사용하며, 안정화가 필요하지 않습니다초수렴 오차 추정 확립 : h k + 1 h^{k+1} h k + 1 후처리 기법 제공 : 초수렴 속도 근사를 구성합니다등가 관계 확립 : 비적합 가상요소법 및 의사응력 형식과의 연결을 보여줍니다유계 영역 Ω ⊂ R d \Omega \subset \mathbb{R}^d Ω ⊂ R d { − Δ u − ∇ p = f in Ω div u = 0 in Ω u = 0 on ∂ Ω \begin{cases}
-\Delta u - \nabla p = f & \text{in } \Omega \\
\text{div } u = 0 & \text{in } \Omega \\
u = 0 & \text{on } \partial\Omega
\end{cases} ⎩ ⎨ ⎧ − Δ u − ∇ p = f div u = 0 u = 0 in Ω in Ω on ∂ Ω
여기서 u u u p p p f f f
H(div)-적합 Raviart-Thomas (RT) 또는 Brezzi-Douglas-Marini (BDM) 요소를 사용합니다:
V ∘ k , ℓ div : = { v h ∈ H 0 ( div , Ω ) : v h ∣ T ∈ P k ( T ; R d ) + H ℓ ( T ) x } \stackrel{\circ}{V}^{\text{div}}_{k,\ell} := \{v_h \in H_0(\text{div},\Omega) : v_h|_T \in P_k(T;\mathbb{R}^d) + H_\ell(T)x\} V ∘ k , ℓ div := { v h ∈ H 0 ( div , Ω ) : v h ∣ T ∈ P k ( T ; R d ) + H ℓ ( T ) x }
여기서 ℓ = k \ell = k ℓ = k ℓ = k − 1 \ell = k-1 ℓ = k − 1
편차 기울기 이산화를 위해 무적 텐서 공간 Σ k t n \Sigma^{tn}_k Σ k t n Σ k t n : = { τ h ∈ Σ k − 1 ( T ) : [ Π F τ n ] F = 0 for all F ∈ F h } \Sigma^{tn}_k := \{\tau_h \in \Sigma^{-1}_k(\mathcal{T}) : [\Pi_F \tau n]_F = 0 \text{ for all } F \in \mathcal{F}_h\} Σ k t n := { τ h ∈ Σ k − 1 ( T ) : [ Π F τ n ] F = 0 for all F ∈ F h }
접선-법선 연속성을 완화하기 위해 다음을 도입합니다:
Λ k = P k ( F h ; R d − 1 ) , Λ ∘ k = P k ( F ∘ h ; R d − 1 ) \Lambda_k = P_k(\mathcal{F}_h;\mathbb{R}^{d-1}), \quad \stackrel{\circ}{\Lambda}_k = P_k(\stackrel{\circ}{\mathcal{F}}_h;\mathbb{R}^{d-1}) Λ k = P k ( F h ; R d − 1 ) , Λ ∘ k = P k ( F ∘ h ; R d − 1 )
약한 편차 기울기 연산자 dev grad w : V ∘ k , ℓ div × Λ ∘ k → Σ k − 1 ( T ) \text{dev grad}_w : \stackrel{\circ}{V}^{\text{div}}_{k,\ell} \times \stackrel{\circ}{\Lambda}_k \to \Sigma^{-1}_k(\mathcal{T}) dev grad w : V ∘ k , ℓ div × Λ ∘ k → Σ k − 1 ( T )
( v , μ ) ∈ H 1 ( T h ; R d ) × L 2 ( F h ; R d − 1 ) (v,\mu) \in H^1(\mathcal{T}_h;\mathbb{R}^d) \times L^2(\mathcal{F}_h;\mathbb{R}^{d-1}) ( v , μ ) ∈ H 1 ( T h ; R d ) × L 2 ( F h ; R d − 1 ) ( dev grad w ( v , μ ) , τ ) T = − ( v , div τ ) T + ( n ⋅ v , n ⊤ τ n ) ∂ T + ( μ , Π F τ n ) ∂ T (\text{dev grad}_w(v,\mu), \tau)_T = -(v, \text{div } \tau)_T + (n \cdot v, n^\top \tau n)_{\partial T} + (\mu, \Pi_F \tau n)_{\partial T} ( dev grad w ( v , μ ) , τ ) T = − ( v , div τ ) T + ( n ⋅ v , n ⊤ τ n ) ∂ T + ( μ , Π F τ n ) ∂ T
u h ∈ V ∘ k , ℓ div u_h \in \stackrel{\circ}{V}^{\text{div}}_{k,\ell} u h ∈ V ∘ k , ℓ div λ h ∈ Λ ∘ k \lambda_h \in \stackrel{\circ}{\Lambda}_k λ h ∈ Λ ∘ k p h ∈ P ℓ ( T h ) / R p_h \in P_\ell(\mathcal{T}_h)/\mathbb{R} p h ∈ P ℓ ( T h ) / R
( dev grad w ( u h , λ h ) , dev grad w ( v h , μ h ) ) + ( div v h , p h ) = ( f , v h ) (\text{dev grad}_w(u_h,\lambda_h), \text{dev grad}_w(v_h,\mu_h)) + (\text{div } v_h, p_h) = (f, v_h) ( dev grad w ( u h , λ h ) , dev grad w ( v h , μ h )) + ( div v h , p h ) = ( f , v h ) ( div u h , q h ) = 0 (\text{div } u_h, q_h) = 0 ( div u h , q h ) = 0
모든 v h ∈ V ∘ k , ℓ div v_h \in \stackrel{\circ}{V}^{\text{div}}_{k,\ell} v h ∈ V ∘ k , ℓ div μ h ∈ Λ ∘ k \mu_h \in \stackrel{\circ}{\Lambda}_k μ h ∈ Λ ∘ k q h ∈ P ℓ ( T h ) / R q_h \in P_\ell(\mathcal{T}_h)/\mathbb{R} q h ∈ P ℓ ( T h ) / R
안정화 없는 설계 : 정교하게 구성된 약한 연산자를 통해 기존 DG 방법에서 필요한 페널티항을 피합니다초수렴성 : 교환 성질을 활용하여 표준 방법보다 한 차수 높은 수렴율을 달성합니다압력 견고성 : 방법이 자연스럽게 압력 견고성을 가집니다저차 가용성 : ( k , ℓ ) = ( 0 , − 1 ) , ( 0 , 0 ) , ( 1 , 0 ) (k,\ell) = (0,-1), (0,0), (1,0) ( k , ℓ ) = ( 0 , − 1 ) , ( 0 , 0 ) , ( 1 , 0 ) 예제 5.1 (2차원) :
정확해: u = curl ψ 2 u = \text{curl } \psi_2 u = curl ψ 2 p = − x 5 − y 5 + 1 3 p = -x^5 - y^5 + \frac{1}{3} p = − x 5 − y 5 + 3 1 여기서 ψ 2 = x 2 ( x − 1 ) 2 y 2 ( y − 1 ) 2 \psi_2 = x^2(x-1)^2y^2(y-1)^2 ψ 2 = x 2 ( x − 1 ) 2 y 2 ( y − 1 ) 2 예제 5.2 (3차원) :
정확해: u = curl ( ψ 3 , ψ 3 , ψ 3 ) T u = \text{curl}(\psi_3, \psi_3, \psi_3)^T u = curl ( ψ 3 , ψ 3 , ψ 3 ) T p = − x 5 − y 5 − z 5 + 1 2 p = -x^5 - y^5 - z^5 + \frac{1}{2} p = − x 5 − y 5 − z 5 + 2 1 여기서 ψ 3 = x 2 ( x − 1 ) 2 y 2 ( y − 1 ) 2 z 2 ( z − 1 ) 2 \psi_3 = x^2(x-1)^2y^2(y-1)^2z^2(z-1)^2 ψ 3 = x 2 ( x − 1 ) 2 y 2 ( y − 1 ) 2 z 2 ( z − 1 ) 2 계산 영역: Ω = ( 0 , 1 ) d \Omega = (0,1)^d Ω = ( 0 , 1 ) d d = 2 , 3 d = 2,3 d = 2 , 3 메시: 균일 심플렉스 메시 분할 구현: MATLAB 패키지 iFEM 기반 속도 오차: ∥ u − u h ∥ \|u - u_h\| ∥ u − u h ∥ 응력 오차: ∥ σ − σ h ∥ \|\sigma - \sigma_h\| ∥ σ − σ h ∥ 압력 오차: ∥ p − p h ∥ \|p - p_h\| ∥ p − p h ∥ 후처리 오차: ∥ u − u h ∗ ∥ \|u - u_h^*\| ∥ u − u h ∗ ∥ ∥ grad h ( u − u h ∗ ) ∥ \|\text{grad}_h(u - u_h^*)\| ∥ grad h ( u − u h ∗ ) ∥ 표 1: 2차원 오차 결과
h (k,ℓ) ∥ u − u h ∥ \|u-u_h\| ∥ u − u h ∥ 차수 ∥ σ − σ h ∥ \|\sigma-\sigma_h\| ∥ σ − σ h ∥ 차수 ∥ p − p h ∥ \|p-p_h\| ∥ p − p h ∥ 차수 2⁻³ (0,0) 2.988e-03 - 3.103e-02 - 7.810e-02 - 2⁻⁴ (0,0) 1.284e-03 1.22 1.677e-02 0.89 3.914e-02 1.00 2⁻⁵ (0,0) 5.988e-04 1.10 8.700e-03 0.95 1.963e-02 1.00
| 2⁻³ | (1,0) | 3.296e-04 | - | 2.447e-03 | - | 7.453e-02 | - |
| 2⁻⁴ | (1,0) | 8.382e-05 | 1.98 | 6.305e-04 | 1.96 | 3.760e-02 | 0.99 |
| 2⁻⁵ | (1,0) | 2.104e-05 | 1.99 | 1.597e-04 | 1.98 | 1.880e-02 | 1.00 |
실험 결과는 이론적 예측을 완전히 검증합니다:
속도 및 응력 : ∥ u − u h ∥ = ∥ σ − σ h ∥ = O ( h k + 1 ) \|u - u_h\| = \|\sigma - \sigma_h\| = O(h^{k+1}) ∥ u − u h ∥ = ∥ σ − σ h ∥ = O ( h k + 1 ) 압력 : ∥ p − p h ∥ = O ( h ) \|p - p_h\| = O(h) ∥ p − p h ∥ = O ( h ) 후처리 속도 : ∥ u − u h ∗ ∥ = O ( h k + 2 ) \|u - u_h^*\| = O(h^{k+2}) ∥ u − u h ∗ ∥ = O ( h k + 2 ) ∥ grad h ( u − u h ∗ ) ∥ = O ( h k + 1 ) \|\text{grad}_h(u - u_h^*)\| = O(h^{k+1}) ∥ grad h ( u − u h ∗ ) ∥ = O ( h k + 1 ) 3차원 실험도 방법의 유효성을 검증하며, 수렴 차수는 이론과 일치합니다.
약한 발산 안정성 조건을 확립합니다:
inf v h ∈ V ∘ k , k − 1 div sup τ h ∈ Σ k t n ( div w τ h , v h ) 0 , h ∥ τ h ∥ div w ∥ v h ∥ = α > 0 \inf_{v_h \in \stackrel{\circ}{V}^{\text{div}}_{k,k-1}} \sup_{\tau_h \in \Sigma^{tn}_k} \frac{(\text{div}_w \tau_h, v_h)_{0,h}}{\|\tau_h\|_{\text{div}_w} \|v_h\|} = \alpha > 0 inf v h ∈ V ∘ k , k − 1 div sup τ h ∈ Σ k t n ∥ τ h ∥ div w ∥ v h ∥ ( div w τ h , v h ) 0 , h = α > 0
정리 : u ∈ H k + 2 ( Ω ; R d ) u \in H^{k+2}(\Omega;\mathbb{R}^d) u ∈ H k + 2 ( Ω ; R d ) ∥ σ − σ h ∥ 0 , h + ∥ dev grad w ( I k , k div u − u h , Q k , F h λ − λ h ) ∥ + ∥ Q ℓ p − p h ∥ ≲ h k + 1 ∣ u ∣ k + 2 \|\sigma - \sigma_h\|_{0,h} + \|\text{dev grad}_w(I^{\text{div}}_{k,k}u - u_h, Q_{k,\mathcal{F}_h}\lambda - \lambda_h)\| + \|Q_\ell p - p_h\| \lesssim h^{k+1}|u|_{k+2} ∥ σ − σ h ∥ 0 , h + ∥ dev grad w ( I k , k div u − u h , Q k , F h λ − λ h ) ∥ + ∥ Q ℓ p − p h ∥ ≲ h k + 1 ∣ u ∣ k + 2
핵심 교환 성질:
Q k t n dev grad = dev grad w I k , k div Q^{tn}_k \text{dev grad} = \text{dev grad}_w I^{\text{div}}_{k,k} Q k t n dev grad = dev grad w I k , k div
이것이 초수렴을 달성하는 핵심입니다.
Scott-Vogelius 요소 : 특수 메시 조건이 필요합니다매끄러운 유한요소 쌍 : 초매끄러운 자유도가 필요합니다분할 메시 위의 적합 쌍 : 구현이 복잡합니다MCS 방법 : k = ℓ ≥ 1 k = \ell \geq 1 k = ℓ ≥ 1 HDG 방법 : 안정화가 필요하지만, 본 논문은 안정화가 필요하지 않습니다가상요소법 : 일반적으로 안정화가 필요하며, 수렴 차수가 낮습니다안정화 없는 초수렴 무발산 유한요소법을 성공적으로 구성했습니다 h k + 1 h^{k+1} h k + 1 방법은 압력 견고성을 가지며 저차 경우를 지원합니다 가상요소법 및 의사응력 형식과의 등가 관계를 확립했습니다 이론적 분석은 주로 다면체 영역을 대상으로 합니다 H² 정칙성 가정이 적용 범위를 제한합니다 구현 복잡도가 상대적으로 높습니다 Navier-Stokes 방정식으로의 확장 적응형 메시 세분화 병렬 알고리즘 개발 공학 응용 검증 이론적 엄밀성 : 완전한 안정성 및 수렴성 분석방법 혁신 : 약한 편차 기울기 연산자 설계가 정교합니다실용적 가치 : 저차 요소를 지원하며 초수렴을 달성합니다충분한 실험 : 2차원 및 3차원 수치 검증이 완전합니다구현 복잡성 : 약한 연산자의 계산이 상대적으로 복잡합니다이론적 제한 : 강한 정칙성 가정이 필요합니다응용 검증 : 실제 공학 문제 검증이 부족합니다본 연구는 Stokes 방정식 수치 방법 분야에서 중요한 이론적 가치를 가지며, 고정밀 무발산 방법 구성을 위한 새로운 아이디어를 제공합니다. 방법의 초수렴 성질과 안정화 없는 특성은 전산유체역학에서 잠재적 응용 가치를 가집니다.
정확한 질량 보존이 필요한 유체 계산 계산 정밀도 요구가 높은 응용 학술 연구 및 방법 검증 논문은 56편의 관련 문헌을 인용하며, 유한요소법, 혼합 방법, 가상요소법 등 다양한 분야의 중요한 연구를 포함하고 있어 문헌 조사가 포괄적입니다.
