혼합유한요소법

편집 역사 토론
수치해석 전산유체역학 구조해석 마지막 수정: 2026-08-30 04:31:52

1. 개요[편집]

혼합유한요소법
Mixed Finite Element Method
정의미지장 두 개 이상을 독립적으로 보간하는 변분 정식화
대표 짝속도–압력 (스토크스) · 플럭스–수두 (다르시) · 변위–응력
선형계안장점 구조 — 대칭 부정부호
안정성 조건커널 위 강제성 + inf-sup(LBB) — 브레치(1974)
위반 증상압력 체커보드 · 해의 비유일성
안정한 예테일러-후드 P2–P1 · MINI · Raviart-Thomas
해법MINRES/GMRES + 블록 전처리 · 우자와 · 혼성화

미지수를 하나 더 늘리면 문제가 어려워지는 게 상식이다. 그런데 가끔은 늘려야 풀린다.

혼합유한요소법(mixed finite element method)은 원래 하나의 미지장으로 환원할 수 있는 문제에서, 그 도함수 성질의 양(응력·플럭스·압력)을 소거하지 않고 독립 미지수로 함께 보간하는 유한요소 정식화다. 포아송 방정식을 예로 들면, 표준 갤러킨 방법은 uu 하나만 두고 q=−∇u\mathbf{q}=-\nabla u 를 사후에 미분으로 뽑는다. 혼합법은 q\mathbf{q} 와 uu 를 동시에 미지수로 두고 1계 연립계

q+∇u=0,∇⋅q=f\mathbf{q} + \nabla u = 0,\qquad \nabla\cdot\mathbf{q} = f

를 각각 약형식으로 만든다. 미지수가 몇 배로 늘고 행렬은 부정부호가 되는데, 그럼에도 이 길을 택하는 이유가 분명히 있다. 그리고 그 대가로 안정성이 공짜가 아니게 된다 — 표준 FEM에서는 요소를 아무거나 골라도 되던 자리에, 두 공간의 궁합을 따지는 inf-sup 조건이 들어선다. 이 문서는 왜 필요한가, 무엇이 부서지는가, 무엇을 고르면 되는가 셋을 다룬다. 안장점 계를 실제로 푸는 반복법 쪽은 우자와 알고리즘이 이미 자세히 다루므로 여기서는 지도만 그린다.

2. 왜 미지수를 늘리는가[편집]

혼합 정식화를 쓰는 동기는 크게 세 가닥이고, 셋이 서로 다른 문제라는 점이 중요하다.

(가) 비압축성 잠금. 선형 탄성에서 포아송비 ν→0.5\nu\to 0.5 로 가면 체적탄성계수가 발산해 ∇⋅u=0\nabla\cdot\mathbf{u}=0 이라는 제약이 사실상 강제된다. 저차 요소는 이 제약을 요소마다 걸면 남는 자유도가 없어져 뻣뻣해진다 — 요소 잠김의 체적 잠김이 이것이다. 정공법은 정수압 pp 를 독립 미지수로 승격시켜 u\mathbf{u}–pp 로 푸는 것이다. 고무 같은 초탄성 재료 모델과 금속 소성(체적 보존 J2J_2 유동)의 상용 요소가 거의 다 혼합형인 이유다. 스토크스 유동은 물리적으로 이 극한이며, 그래서 비압축 탄성과 스토크스는 같은 수학이다.

(나) 플럭스의 국소 보존. 다르시 법칙 문서가 지적하듯, 수두를 조각별 다항식으로 잡으면 q=−K∇h\mathbf{q}=-\mathbf{K}\nabla h 는 요소 경계에서 법선 성분이 불연속이다. 셀에 들어간 물과 나간 물이 안 맞고, 그 차이가 뒤따르는 대류-확산 방정식 솔버에게는 없는 우물과 없는 샘으로 보인다. 오염물 농도가 음수로 튀는 고전적 사고가 여기서 난다. 혼합법은 q\mathbf{q} 를 법선 연속성이 자유도 정의에 박힌 공간에서 뽑아 이 문제를 구조적으로 없앤다.

(다) 도함수 양의 정확도. 응력·모멘트가 최종 산출물인데 변위만 보간하면, 응력은 한 차수 낮은 정확도로 나오고 요소 경계에서 불연속이다. 응력을 직접 보간하면 그 차수를 되찾는다. 판·셸의 라이스너-미틀린 정식화에서 전단력을 독립 보간해 전단 잠김을 피하는 것도 이 갈래다.

3. 안장점 구조[편집]

두 미지장을 이산화하면 선형계가 언제나 다음 블록 꼴이 된다.

(ABTB0)(up)=(fg)\begin{pmatrix} A & B^{\mathsf T}\\ B & 0\end{pmatrix} \begin{pmatrix}\mathbf{u}\\ \mathbf{p}\end{pmatrix} = \begin{pmatrix}\mathbf{f}\\ \mathbf{g}\end{pmatrix}

AA 는 대칭 양반정부호(점성·강성·투수 저항), BB 는 제약 연산자(이산 발산), p\mathbf{p} 는 그 제약의 승수다. 이 계의 성격을 규정하는 사실 하나만 기억하면 된다. 대각에 0 블록이 있으므로 이 행렬은 절대 양정부호가 아니다. AA 가 ker⁡B\ker B 위에서 양정부호이고 BB 가 행 계수 full 이면, 관성(고윳값 부호 개수)은 정확히 속도 자유도 수만큼의 양의 고윳값과 압력 자유도 수만큼의 음의 고윳값이다.

실무적 귀결은 즉각적이다.

  • 켤레기울기법을 그대로 못 쓴다. 대칭이지만 부정부호이므로 MINRES 또는 GMRES 계열로 가거나, 압력을 소거한 슈어 보수 S=BA−1BTS = BA^{-1}B^{\mathsf T} 로 내려가야 한다. 후자가 우자와 알고리즘의 세계다.
  • 직접법도 조심해야 한다. 촐레스키 분해는 아예 안 되고, LDLTLDL^{\mathsf T} 를 쓰더라도 2×22\times2 피봇을 허용하는 대칭 부정부호 전용 인수분해(Bunch-Kaufman 계열)가 필요하다. 대각에 0이 들어오는 순간 단순 피봇팅이 무너지기 때문이다.

4. LBB(inf-sup) 조건 — 궁합의 문제[편집]

혼합법이 표준 FEM과 근본적으로 다른 지점이 여기다. 표준 갤러킨에서는 공간을 촘촘하게만 하면 수렴이 따라오지만, 혼합법에서는 두 공간의 짝이 맞아야 한다. 브레치(F. Brezzi, 1974)가 정리한 두 조건은 다음과 같다.

  1. 커널 위 강제성. a(⋅,⋅)a(\cdot,\cdot) 가 ker⁡Bh={vh:b(vh,qh)=0 ∀qh}\ker B_h = \{v_h : b(v_h,q_h)=0\ \forall q_h\} 위에서 hh 에 무관한 상수로 강제적(coercive)일 것.
  2. 이산 inf-sup 조건. 격자 크기 hh 에 무관한 β>0\beta>0 가 존재해
inf⁡qh≠0 sup⁡vh≠0 b(vh,qh)∥vh∥V ∥qh∥Q ≥ β>0\inf_{q_h\neq 0}\ \sup_{v_h\neq 0}\ \frac{b(v_h,q_h)}{\lVert v_h\rVert_V\,\lVert q_h\rVert_Q}\ \ge\ \beta > 0

읽는 법은 이렇다. 어떤 압력 모드를 가져와도, 그것을 “느끼는” 속도 함수가 이산 공간 안에 충분히 크게 존재해야 한다. 존재하지 않는 압력 모드 qhq_h 가 있다면 BTqh=0B^{\mathsf T}q_h = 0 이므로 그 모드는 해에 아무 흔적을 남기지 않은 채 자유롭게 더해질 수 있고, 압력은 유일하지 않게 된다.

여기서 가장 흔한 오해를 짚어야 한다. 연속 문제의 inf-sup 조건이 성립한다는 사실은 이산 조건을 전혀 보장하지 않는다. 연속 공간에서는 발산 연산자가 전사(surjective)라 조건이 성립하지만(H01\mathrm{H}^1_0 에서 평균 0인 L2L^2 로), 유한차원으로 자르면 속도 공간이 압력 공간에 비해 “좁아져” 전사성이 깨질 수 있다. 이산화는 부분공간을 취하는 일이고, 부분공간을 취하면 sup 가 줄어든다. 혼합법이 어려운 이유의 절반이 이 한 줄이다.

그리고 β\beta 는 존재만으로는 부족하고 hh 에 무관해야 한다. βh→0\beta_h \to 0 이면 오차 상한 ∥p−ph∥≤C(1+1/βh)inf⁡∥p−qh∥\lVert p - p_h\rVert \le C(1+1/\beta_h)\inf\lVert p-q_h\rVert 가 격자를 조일수록 나빠지고, 우자와 알고리즘 문서가 보여주듯 슈어 보수의 조건수 Γ2/β2\Gamma^2/\beta^2 가 폭발해 반복 솔버까지 같이 죽는다. 이산화 이론의 상수가 솔버 성능으로 직접 번역되는 드문 사례다.

5. 압력 체커보드 — 위반의 얼굴[편집]

LBB를 위반하면 압력장에 인접 절점마다 부호가 뒤집히는 체커보드(checkerboard) 진동이 나타난다. 정사각 격자 위의 Q1\mathrm{Q}_1 속도 – P0\mathrm{P}_0 압력 쌍이 교과서적 사례로, 셀마다 (−1)i+j(-1)^{i+j} 로 값을 주는 압력 모드가 이산 발산 연산자의 커널에 정확히 들어간다. 속도해는 멀쩡한데 압력만 못 쓰게 되는 것이 특징이며, 등차수 쌍(Q1\mathrm{Q}_1–Q1\mathrm{Q}_1, P1\mathrm{P}_1–P1\mathrm{P}_1)에서는 아예 여러 개의 가짜 모드가 뜬다.

체적 잠김에서 정수압이 요소마다 진동하는 현상과 같은 뿌리라는 점이 재미있다. 자유도 대비 제약이 너무 많으면 잠기고(속도가 죽고), 너무 적으면 체커보드가 뜬다(압력이 죽는다). LBB는 이 줄다리기의 정확한 균형점을 진술한 것이다.

경계할 관행이 하나 있다. 체커보드가 보인다고 압력을 후처리로 평활화(smoothing)하는 것. 그림은 예뻐지지만 불안정한 이산화가 만든 오차가 사라지는 것이 아니라 눈에 안 보이게 될 뿐이다. 진단으로 쓰는 것은 좋고 처방으로 쓰면 안 된다.1

자유도 개수를 세어 보는 constraint counting 휴리스틱(속도 자유도 대 압력 자유도의 비)이 널리 쓰이지만, 이건 필요조건 감각일 뿐 충분조건이 아니다. Q1\mathrm{Q}_1–P0\mathrm{P}_0 는 개수만 보면 통과하는데도 체커보드 모드가 뜬다. 최종 판정은 이론 증명이거나 아래의 수치 검사다.

6. 안정한 요소쌍[편집]

짝무대성격
테일러-후드 P2–P1 (Q2–Q1)스토크스, 비압축 탄성연속 압력. 가장 널리 쓰이는 국룰
MINI (P1+버블 – P1)스토크스속도에 버블 하나 얹어 저차 유지
크루제-라비아르 (P2+버블 – P1 불연속)스토크스불연속 압력, 요소별 질량 보존
Q2–P1 불연속CFD 상용 코드국소 보존이 필요할 때
Raviart-Thomas RTk다르시, 확산H(div) 적합, 플럭스가 주역
BDMk다르시RT보다 완전한 다항식, 같은 법선 차수

테일러-후드 요소(1973)의 처방은 단순하다 — 속도를 압력보다 한 차수 높게 잡는다. 압력이 연속이라 자유도가 싸고, 대부분의 격자에서 β\beta 가 hh 에 무관함이 증명되어 있다. 단점은 압력이 연속이라 요소별 질량 보존이 정확히는 성립하지 않는다는 것.

MINI 요소(Arnold-Brezzi-Fortin, 1984)는 반대 전략이다. P1\mathrm{P}_1–P1\mathrm{P}_1 은 불안정하지만, 속도에 요소 내부에서만 살고 경계에서 0인 버블 함수를 하나 더해 주면 sup 가 딱 필요한 만큼 커져 조건을 통과한다. 버블은 정적 축약으로 요소 단계에서 지워지므로 전역 자유도가 거의 안 는다. “안정성만 사고 비용은 안 낸다”는 점에서 우아한 트릭이며, 실제로 버블 축약 후의 계가 뒤에 나올 안정화 항과 형태적으로 유사해진다는 것이 알려져 있다.

Raviart-Thomas(1977, 3차원 확장은 네델렉)와 BDM(Brezzi-Douglas-Marini, 1985)은 다르시 계열의 표준이다. 핵심은 이들이 H(div)H(\mathrm{div}) 적합 공간이라는 것 — 자유도가 면 위의 법선 성분 적분으로 정의되므로 면을 가로지르는 법선 플럭스가 자동으로 연속이고, 발산이 압력 공간으로 정확히 사상되어 ∇⋅qh=fh\nabla\cdot\mathbf{q}_h = f_h 가 요소마다 정확히 성립한다. 국소 질량 보존이 근사가 아니라 항등식이 된다는 뜻이며, 이송 솔버에 넘길 속도장으로는 이만한 것이 없다. 최저차 RT0\mathrm{RT}_0 를 사각격자에서 쓰면 압력이 셀 중심에서 초수렴하고, 혼성화 후의 계가 유한체적법의 2점 플럭스 도식과 놀랍도록 닮는다.2

RT와 BDM의 차이는 한 줄로 정리된다. RTk\mathrm{RT}_k 는 완전 다항식 [Pk]d[\mathbb{P}_k]^d 에 특정 방향 항을 조금 더한 공간이고, BDMk\mathrm{BDM}_k 는 [Pk]d[\mathbb{P}_k]^d 전체를 쓴다. 같은 법선 차수에서 BDM 쪽이 플럭스 근사 차수가 높지만 발산의 상은 한 차수 낮다. 즉 플럭스 정확도를 살 것인가, 발산 쪽 균형을 살 것인가의 거래다.

7. 안정화 — 등차수 쌍을 살리는 길[편집]

LBB를 만족하는 쌍만 쓰기엔 불편한 사정이 많다. 속도와 압력의 자유도 배치가 다르면 자료구조가 복잡해지고, 코드 재사용도 어렵다. 그래서 불안정한 등차수 쌍을 쓰되 식을 고쳐서 살리는 노선이 발달했다.

  • PSPG(Pressure-Stabilizing Petrov-Galerkin, Hughes-Franca-Balestra 1986). 논문 제목 자체가 “바부슈카-브레치 조건을 우회하기”다. 압력 시험함수에 τ∇qh\tau\nabla q_h 항을 더해 잔차에 가중을 실으면, 블록계의 (2,2)(2,2) 자리에 −C-C 가 생겨 안장점이 아니라 정칙화된 안장점이 된다. 이 CC 가 사실상 이산 압력 라플라시안이라, 체커보드처럼 고주파로 진동하는 압력 모드를 정확히 겨냥해 감쇠시킨다.
  • GLS(Galerkin/Least-Squares, Hughes-Franca-Hulbert 1989). SUPG와 PSPG를 잔차 최소제곱이라는 하나의 틀로 통합한 것. 대류 안정화와 압력 안정화를 한 매개변수 체계로 다룬다.
  • 국소 사영 안정화·압력 사영법(Bochev-Dohrmann-Gunzburger 2006 계열). 잔차 대신 “압력에서 요소별 평균을 뺀 성분”에만 벌점을 주는 방식. 잔차 기반이 아니어서 시간항·비선형항의 미분을 안 챙겨도 되고, 매개변수 민감도가 낮다는 것이 장점이다.

안정화의 대가는 매개변수 τ\tau 다. 너무 작으면 체커보드가 안 잡히고, 너무 크면 해가 과도하게 뭉개지며, 적정 τ\tau 는 격자 크기·물성·시간간격에 의존한다. “안정한 요소쌍을 쓰면 튜닝 노브가 없고, 안정화를 쓰면 노브가 생긴다”는 것이 이 선택의 본질이다. 상용 CFD 코드가 등차수 + PSPG 조합을 선호하는 이유는 성능이 아니라 구현 단순성 쪽에 가깝다.

8. 푸는 법 — 지도만[편집]

블록계를 실제로 푸는 방법은 우자와 알고리즘·증강 라그랑주법·크리로프 부분공간법 문서에 흩어져 있다. 여기서는 좌표만 찍는다.

  • 슈어 보수 경로. u\mathbf{u} 를 소거해 Sp=g^S\mathbf{p}=\hat{\mathbf{g}} 로 내려간다. SS 는 조립하지 않고 곱셈만 쓴다. 리처드슨으로 밀면 우자와, CG로 밀면 슈어 보수 CG. 전처리기는 압력 질량행렬 MpM_p 가 정석인데, 그 근거가 위의 β2≤(Sq,q)/(Mpq,q)≤Γ2\beta^2 \le (Sq,q)/(M_pq,q) \le \Gamma^2 부등식이다. 이산화 상수가 전처리기를 지목해 주는 셈.
  • 블록 전처리 경로. 소거하지 않고 전체 계를 MINRES/GMRES에 넣되, 전처리기의 (1,1)(1,1) 블록에 AA 의 다중격자법 근사를, (2,2)(2,2) 블록에 SS 의 근사를 꽂는다. 대규모 스토크스·나비에-스토크스 솔버의 현대 표준.
  • 혼성화 경로. RT/BDM 같은 H(div)H(\mathrm{div}) 쌍에서는 면 법선 연속성을 자유도로 강제하는 대신 면마다 라그랑주 승수를 하나 더 두고 요소 내부 자유도를 정적 축약한다. 남는 것은 면 승수만의 계이고, 이게 대칭 양정부호다. 부정부호 문제가 정부호 문제로 바뀌므로 CG와 다중격자를 그대로 태울 수 있다. 미지수를 더 늘려서 문제를 쉽게 만든다는, 이 문서 전체의 주제를 한 번 더 반복하는 트릭이며, 불연속 갤러킨법 계열의 HDG로 이어진다.

9. 새 요소쌍을 만들었다면 — 수치 inf-sup 검사[편집]

이론 증명 전에 실무자가 돌리는 표준 절차가 있다. 일반화 고유값 문제

BMu−1BT q=λ Mp qB M_u^{-1} B^{\mathsf T}\,\mathbf{q} = \lambda\, M_p\,\mathbf{q}

를 풀어 0이 아닌 최소 고유값을 찾고 βh=λmin⁡\beta_h=\sqrt{\lambda_{\min}} 로 둔다. 격자를 단계적으로 조이면서 βh\beta_h 가 상수로 수렴하면 안정, 0으로 떨어지면 불안정이다. 정확히 0인 고유값의 개수가 가짜 압력 모드의 개수이며, 그 고유벡터를 그려 보면 체커보드가 눈으로 보인다. 실습 과제로도 좋은 검사다.

같은 판정을 우자와 반복 수로 대신할 수도 있다 — 격자를 절반으로 줄일 때마다 반복 수가 뚜렷하게 늘면 βh→0\beta_h\to 0 이라는 강력한 정황 증거다.

10. 여담 — 잠금 대책은 사실 다 혼합법이었다[편집]

요소 잠김 문서가 나열하는 처방들 — 선택적 감소적분(SRI), B-bar — 은 겉보기엔 “적분점을 덜 쓰는 요령”이다. 그런데 말커스와 휴즈(1978)가 보인 결과가 유명하다. 적절한 조건에서 선택적 감소적분은 불연속 압력을 쓰는 혼합 정식화와 수학적으로 동치다. 즉 현장에서 “적분점 하나 뺐더니 잠김이 풀렸다”고 말할 때, 실제로 일어난 일은 압력을 요소별 상수로 독립 보간한 혼합법으로 갈아탄 것이었다.

그래서 감소적분 요소가 왜 어떤 격자에서 갑자기 체커보드 압력을 뱉는지도 같은 언어로 설명된다 — 그 등가 혼합 쌍이 LBB를 만족하지 않기 때문이다. 요령처럼 보이던 것에 이름과 이론이 붙는 순간, 언제 깨지는지도 같이 알게 된다. 수치해석에서 “왜 되는지 알아내는” 작업의 값어치가 대개 여기에 있다.3

11. 관련 문서[편집]

12. Footnotes[편집]

  1. 압력 평활화가 무해한 경우도 있긴 하다 — 이미 안정한 쌍을 썼는데 초기 반복의 잔여 진동을 보기 좋게 만드는 정도. 문제는 보고서를 받는 쪽이 그 차이를 모른다는 것이다. 체커보드는 그림의 문제가 아니라 이산화가 답을 유일하게 정하지 못하고 있다는 신호이므로, 지우기 전에 요소쌍부터 확인하는 것이 순서다. ↩

  2. 그래서 RT0\mathrm{RT}_0 혼합법과 셀 중심 유한체적을 “사실상 같은 도식”이라 부르는 문헌이 많다. 정확히는 사각격자·대각 이방성 같은 좋은 조건에서 질량행렬을 대각화(mass lumping)했을 때 일치한다. 격자가 뒤틀리거나 K\mathbf{K} 의 비주축 성분이 커지면 갈라지고, 이때 정합성을 지키는 쪽은 혼합법이다. ↩

  3. 반대 방향의 교훈도 있다. 이론 없이 요령만 물려받은 코드는 조건이 바뀌는 순간 조용히 틀린다. 감소적분 요소를 준 대로 쓰다가 압력이 이상해졌을 때, “적분점을 뺐다”는 설명으로는 원인에 도달할 수 없지만 “등가 쌍이 LBB를 어긴다”는 설명으로는 도달한다. 같은 현상에 대한 두 서술의 값이 이만큼 다르다. ↩