6.155 전처리 기법과 수렴 가속
1. 전처리의 기본 원리
반복법의 수렴 속도는 계수 행렬 \mathbf{A}의 조건수 \kappa(\mathbf{A})에 크게 의존한다. 켤레 기울기법(CG)의 경우, 오차를 \epsilon 이하로 줄이는 데 필요한 반복 횟수는 \mathcal{O}(\sqrt{\kappa(\mathbf{A})})에 비례한다. 따라서 조건수가 큰 시스템에서는 수렴이 매우 느려질 수 있다.
전처리(preconditioning)는 원래 시스템 \mathbf{A}\mathbf{x} = \mathbf{b}를 등가이면서 조건수가 더 작은 시스템으로 변환하는 기법이다. 전처리 행렬(preconditioner) \mathbf{P} \approx \mathbf{A}를 도입하면, 전처리된 시스템은 다음과 같다.
\mathbf{P}^{-1}\mathbf{A}\mathbf{x} = \mathbf{P}^{-1}\mathbf{b} \quad (\text{좌전처리, left preconditioning})
\mathbf{A}\mathbf{P}^{-1}\mathbf{u} = \mathbf{b}, \quad \mathbf{x} = \mathbf{P}^{-1}\mathbf{u} \quad (\text{우전처리, right preconditioning})
\mathbf{P}_L^{-1}\mathbf{A}\mathbf{P}_R^{-1}\mathbf{u} = \mathbf{P}_L^{-1}\mathbf{b}, \quad \mathbf{x} = \mathbf{P}_R^{-1}\mathbf{u} \quad (\text{양측 전처리, split preconditioning})
이상적인 전처리 행렬은 다음 조건을 만족하여야 한다.
- \kappa(\mathbf{P}^{-1}\mathbf{A}) \ll \kappa(\mathbf{A}): 전처리된 시스템의 조건수가 원래보다 훨씬 작아야 한다.
- \mathbf{P}\mathbf{z} = \mathbf{r}의 풀이 비용이 낮아야 한다: 매 반복마다 전처리 시스템을 풀어야 하므로, 이 비용이 크면 전처리의 이점이 상쇄된다.
- \mathbf{P}의 구성 비용이 합리적이어야 한다.
전처리된 켤레 기울기법(PCG)
대칭 양의 정부호(SPD) 행렬 \mathbf{A}에 대한 전처리된 CG(Preconditioned Conjugate Gradient, PCG)의 알고리즘은 다음과 같다.
초기화: \mathbf{r}_0 = \mathbf{b} - \mathbf{A}\mathbf{x}_0, \mathbf{z}_0 = \mathbf{P}^{-1}\mathbf{r}_0, \mathbf{p}_0 = \mathbf{z}_0
k = 0, 1, 2, \ldots에 대해 반복한다.
\alpha_k = \frac{\mathbf{r}_k^\top \mathbf{z}_k}{\mathbf{p}_k^\top \mathbf{A}\mathbf{p}_k}
\mathbf{x}_{k+1} = \mathbf{x}_k + \alpha_k \mathbf{p}_k
\mathbf{r}_{k+1} = \mathbf{r}_k - \alpha_k \mathbf{A}\mathbf{p}_k
\mathbf{z}_{k+1} = \mathbf{P}^{-1}\mathbf{r}_{k+1}
\beta_k = \frac{\mathbf{r}_{k+1}^\top \mathbf{z}_{k+1}}{\mathbf{r}_k^\top \mathbf{z}_k}
\mathbf{p}_{k+1} = \mathbf{z}_{k+1} + \beta_k \mathbf{p}_k
수렴 판정: \|\mathbf{r}_{k+1}\| / \|\mathbf{b}\| < \epsilon이면 종료한다.
PCG의 수렴 상한은 전처리된 시스템의 조건수에 의해 결정된다.
\frac{\|\mathbf{x}_k - \mathbf{x}^*\|_\mathbf{A}}{\|\mathbf{x}_0 - \mathbf{x}^*\|_\mathbf{A}} \leq 2\left(\frac{\sqrt{\kappa(\mathbf{P}^{-1}\mathbf{A})} - 1}{\sqrt{\kappa(\mathbf{P}^{-1}\mathbf{A})} + 1}\right)^k
2. 대각 전처리와 야코비 전처리
가장 간단한 전처리는 대각 전처리(diagonal preconditioning) 또는 야코비 전처리(Jacobi preconditioning)이다.
\mathbf{P} = \text{diag}(a_{11}, a_{22}, \ldots, a_{nn})
이 전처리의 장점은 구성과 적용이 \mathcal{O}(n)이라는 점이다. \mathbf{P}^{-1}\mathbf{r}의 계산은 성분별 나눗셈으로 수행된다.
z_i = r_i / a_{ii}, \quad i = 1, 2, \ldots, n
대각 전처리는 행렬의 대각 원소 크기가 크게 다른 경우(스케일링 불균형)에 효과적이다. 로봇공학에서 관성 행렬 \mathbf{M}(\mathbf{q})은 관절의 종류(회전 관절 vs. 병진 관절)와 링크 질량에 따라 대각 원소의 크기가 수십 배 이상 차이날 수 있으므로, 대각 전처리만으로도 유의미한 수렴 개선을 얻을 수 있다.
3. 불완전 분해 전처리
3.1 불완전 LU 분해(ILU)
불완전 LU 분해(Incomplete LU, ILU)는 LU 분해 과정에서 특정 위치의 채움(fill-in)을 허용하지 않아 희소성을 유지하는 방법이다.
ILU(0)은 원래 행렬 \mathbf{A}의 희소 패턴과 동일한 패턴만을 유지한다. 즉, \mathbf{A}에서 비영인 위치에서만 \mathbf{L}과 \mathbf{U}의 원소를 허용한다.
\mathbf{A} \approx \tilde{\mathbf{L}}\tilde{\mathbf{U}}, \quad \text{단, } (\tilde{\mathbf{L}} + \tilde{\mathbf{U}})_{ij} = 0 \text{ if } a_{ij} = 0
ILU(p)는 레벨-p 채움을 허용하여 더 정확한 근사를 제공한다. p가 증가할수록 전처리 품질은 향상되지만 비용도 증가한다.
ILU 전처리의 적용은 삼각 풀이(triangular solve)를 통해 수행된다.
\mathbf{P}\mathbf{z} = \mathbf{r} \quad \Longrightarrow \quad \tilde{\mathbf{L}}\tilde{\mathbf{U}}\mathbf{z} = \mathbf{r}
이를 전방 대입과 후방 대입으로 \mathcal{O}(\text{nnz})에 풀 수 있다.
3.2 불완전 Cholesky 분해(IC)
SPD 행렬에 대해서는 불완전 Cholesky 분해(Incomplete Cholesky, IC)를 사용할 수 있다.
\mathbf{A} \approx \tilde{\mathbf{L}}\tilde{\mathbf{L}}^\top
IC(0)은 ILU(0)의 대칭 버전으로, 메모리 사용량이 ILU(0)의 약 절반이며, SPD 성질을 보존하므로 PCG와의 호환성이 보장된다. IC 분해의 각 원소는 다음과 같이 계산된다.
\tilde{l}_{ii} = \sqrt{a_{ii} - \sum_{k=1}^{i-1} \tilde{l}_{ik}^2}
\tilde{l}_{ij} = \frac{1}{\tilde{l}_{jj}}\left(a_{ij} - \sum_{k=1}^{j-1} \tilde{l}_{ik}\tilde{l}_{jk}\right), \quad i > j
단, a_{ij} = 0인 위치에서는 \tilde{l}_{ij} = 0으로 설정한다.
SSOR 전처리
대칭 연속 과완화법(Symmetric Successive Over-Relaxation, SSOR)을 전처리로 활용할 수 있다. 행렬 \mathbf{A} = \mathbf{D} + \mathbf{L} + \mathbf{L}^\top (SPD 행렬의 경우)에 대해 SSOR 전처리 행렬은 다음과 같다.
\mathbf{P}_{\text{SSOR}} = \frac{1}{\omega(2 - \omega)}(\mathbf{D} + \omega\mathbf{L})\mathbf{D}^{-1}(\mathbf{D} + \omega\mathbf{L}^\top)
\omega = 1이면 대칭 가우스-자이델 전처리가 된다. SSOR 전처리의 장점은 별도의 전처리 구성 단계가 불필요하며, 행렬 \mathbf{A}의 원소만으로 직접 적용할 수 있다는 점이다.
4. 다중 격자 전처리
4.1 기본 원리
다중 격자(multigrid) 방법은 여러 해상도의 격자(grid)를 사용하여 반복법의 수렴을 가속시키는 기법이다. 고전적 반복법(야코비, 가우스-자이델)은 고주파 오차 성분은 빠르게 감쇄시키지만 저주파 오차 성분은 느리게 감쇄시킨다. 다중 격자 방법은 이 저주파 오차를 거친 격자(coarse grid)에서 효율적으로 제거한다.
4.2 V-사이클 다중 격자
가장 기본적인 다중 격자 구조인 V-사이클은 다음 단계로 구성된다.
- 전 스무딩(pre-smoothing): 미세 격자(fine grid)에서 \nu_1번의 스무딩 반복(예: 가우스-자이델)을 수행한다.
- 제한(restriction): 잔차를 거친 격자로 전이한다. \mathbf{r}^H = \mathbf{I}_h^H \mathbf{r}^h
- 거친 격자 풀이: 거친 격자에서 오차 방정식 \mathbf{A}^H \mathbf{e}^H = \mathbf{r}^H를 풀이한다. (재귀적으로 더 거친 격자로 진행할 수 있다.)
- 보간(prolongation): 거친 격자의 보정값을 미세 격자로 전이한다. \mathbf{e}^h = \mathbf{I}_H^h \mathbf{e}^H
- 보정(correction): \mathbf{x}^h \leftarrow \mathbf{x}^h + \mathbf{e}^h
- 후 스무딩(post-smoothing): 미세 격자에서 \nu_2번의 스무딩 반복을 수행한다.
여기서 \mathbf{I}_h^H는 제한 연산자, \mathbf{I}_H^h는 보간 연산자이다.
4.3 대수적 다중 격자(AMG)
기하학적 다중 격자는 격자의 기하학적 구조를 알아야 하지만, 대수적 다중 격자(Algebraic Multigrid, AMG)는 행렬 \mathbf{A}의 원소값만으로 거친 격자를 자동으로 구성한다. AMG의 핵심 단계는 다음과 같다.
- 강한 연결(strong connection) 판정: \vert a_{ij} \vert \geq \theta \max_{k \neq i} \vert a_{ik} \vert이면 i와 j가 강하게 연결되었다고 판정한다. 여기서 \theta \in (0, 1)은 강도 임계값이다.
- 거친 격자 선택(coarsening): 강한 연결 관계를 이용하여 거친 격자 점을 선택한다.
- 보간 연산자 구성: 미세 격자와 거친 격자 간의 보간 가중치를 결정한다.
AMG를 PCG의 전처리로 사용하면 \mathcal{O}(n)에 가까운 총 계산량으로 수렴을 달성할 수 있어, 대규모 희소 SPD 시스템에 대해 최적에 가까운 복잡도를 보인다.
5. 블록 전처리
로봇공학에서 나타나는 구조화된 선형 시스템에는 블록 전처리(block preconditioning)가 효과적이다.
5.1 블록 대각 전처리
다관절 로봇에서 관절 간의 결합이 약한 경우, 관성 행렬의 블록 대각 부분을 전처리로 사용할 수 있다.
\mathbf{P}_{\text{block}} = \text{blkdiag}(\mathbf{A}_{11}, \mathbf{A}_{22}, \ldots, \mathbf{A}_{kk})
각 블록 \mathbf{A}_{ii}는 소규모이므로 직접법(Cholesky 분해)으로 효율적으로 분해할 수 있으며, 전처리 적용은 각 블록에 대한 삼각 풀이로 수행된다.
슈어 보수 전처리
새들 포인트(saddle point) 구조를 가지는 시스템에서는 슈어 보수(Schur complement) 기반 전처리가 효과적이다. 다음과 같은 블록 시스템을 고려한다.
\begin{bmatrix} \mathbf{A} & \mathbf{B}^\top \\ \mathbf{B} & \mathbf{0} \end{bmatrix} \begin{bmatrix} \mathbf{x} \\ \boldsymbol{\lambda} \end{bmatrix} = \begin{bmatrix} \mathbf{f} \\ \mathbf{g} \end{bmatrix}
이 시스템은 구속 조건이 있는 로봇 동역학이나 접촉 문제에서 빈번히 나타난다. 슈어 보수 \mathbf{S} = \mathbf{B}\mathbf{A}^{-1}\mathbf{B}^\top에 대한 근사 전처리를 구성하면 효율적인 블록 전처리를 얻을 수 있다.
블록 삼각 전처리의 한 형태는 다음과 같다.
\mathbf{P} = \begin{bmatrix} \mathbf{A} & \mathbf{0} \\ \mathbf{B} & -\tilde{\mathbf{S}} \end{bmatrix}
여기서 \tilde{\mathbf{S}} \approx \mathbf{S}는 슈어 보수의 근사이다. \tilde{\mathbf{S}}가 \mathbf{S}를 잘 근사할수록 전처리 품질이 향상된다.
전처리 기법의 비교
| 전처리 방법 | 구성 비용 | 적용 비용 | 메모리 | 효과 |
|---|---|---|---|---|
| 야코비(대각) | \mathcal{O}(n) | \mathcal{O}(n) | \mathcal{O}(n) | 약함 |
| SSOR | \mathcal{O}(1) | \mathcal{O}(\text{nnz}) | \mathcal{O}(1) | 중간 |
| IC(0) | \mathcal{O}(\text{nnz}) | \mathcal{O}(\text{nnz}) | \mathcal{O}(\text{nnz}) | 중간~양호 |
| ILU(0) | \mathcal{O}(\text{nnz}) | \mathcal{O}(\text{nnz}) | \mathcal{O}(\text{nnz}) | 중간~양호 |
| ILU(p), p > 0 | > \mathcal{O}(\text{nnz}) | > \mathcal{O}(\text{nnz}) | > \mathcal{O}(\text{nnz}) | 양호 |
| AMG | \mathcal{O}(n) ~ \mathcal{O}(n \log n) | \mathcal{O}(n) | \mathcal{O}(\text{nnz}) | 우수 |
수렴 가속의 이론적 분석
전처리가 수렴에 미치는 효과를 정량적으로 분석한다. 원래 시스템의 조건수가 \kappa(\mathbf{A})이고 전처리 후 \kappa(\mathbf{P}^{-1}\mathbf{A})가 되었다고 하면, PCG에서 잔차를 \epsilon 이하로 줄이기 위한 반복 횟수의 비는 대략 다음과 같다.
\frac{k_{\text{PCG}}}{k_{\text{CG}}} \approx \sqrt{\frac{\kappa(\mathbf{P}^{-1}\mathbf{A})}{\kappa(\mathbf{A})}}
예를 들어, \kappa(\mathbf{A}) = 10^6이고 IC(0) 전처리로 \kappa(\mathbf{P}^{-1}\mathbf{A}) = 10^2로 줄였다면, 필요 반복 횟수는 약 \sqrt{10^2 / 10^6} = 1/100 배로 감소한다. 즉, CG에서 1000번 반복이 필요했던 문제를 PCG에서 10번 반복으로 풀 수 있게 된다.
6. 로봇공학에서의 전처리 적용 사례
6.1 다체 동역학 시뮬레이션
다수의 강체(rigid body)가 접촉하는 시뮬레이션에서, 접촉력 계산을 위한 선형 시스템의 조건수는 접촉 강성(contact stiffness)에 따라 10^6 이상으로 커질 수 있다. 대각 스케일링과 불완전 분해를 결합한 전처리가 실시간 시뮬레이션에서 널리 사용된다.
6.2 유한 요소 해석
소프트 로봇의 변형 해석에서 강성 행렬 \mathbf{K}는 대규모 희소 SPD 행렬이며, 재료 특성의 불균일성으로 인해 조건수가 매우 클 수 있다. AMG 전처리된 PCG는 이러한 시스템에서 격자 크기에 거의 독립적인 반복 횟수를 보이므로, 대규모 유한 요소 문제의 표준 풀이 방법으로 자리 잡았다.
6.3 궤적 최적화
전신 인간형 로봇의 궤적 최적화에서 나타나는 KKT 시스템은 시간 단계에 걸친 블록 삼대각(block tridiagonal) 구조를 가진다. 이 구조를 이용한 블록 전처리는 문제의 시간 지평(time horizon)에 대해 선형적인 계산량으로 풀이를 가능하게 한다.
참고문헌
- Saad, Y. (2003). Iterative Methods for Sparse Linear Systems (2nd ed.). SIAM.
- Golub, G. H., & Van Loan, C. F. (2013). Matrix Computations (4th ed.). Johns Hopkins University Press.
- Benzi, M., Golub, G. H., & Liesen, J. (2005). “Numerical solution of saddle point problems.” Acta Numerica, 14, 1–137.
- Briggs, W. L., Henson, V. E., & McCormick, S. F. (2000). A Multigrid Tutorial (2nd ed.). SIAM.
- Trefethen, L. N., & Bau, D. (1997). Numerical Linear Algebra. SIAM.
v 0.1