실베스터 방정식

편집 역사 토론
수치해석 최적설계 마지막 수정: 2026-08-10 04:33:18

1. 개요[편집]

실베스터 방정식
Sylvester equation
형태$AX + XB = C$ ($X$가 미지 행렬)
크기$A: m\times m$, $B: n\times n$, $X, C: m \times n$
유일해 조건$\sigma(A) \cap \sigma(-B) = \varnothing$
제안J. J. Sylvester, 1884
표준 해법바텔스-스튜어트 알고리즘, $O(m^3+n^3)$
특수 경우리아푸노프 방정식 ($B = A^{\top}$)

미지수가 스칼라면 방정식, 벡터면 선형계, 행렬이면 — 그것도 그냥 선형계다. 다만 눈에 안 보이게 크다.

실베스터 방정식(Sylvester equation)은 미지 행렬 XX에 대한 선형 행렬방정식

AX+XB=CAX + XB = C

이다. AAm×mm\times m, BBn×nn\times n, XXCCm×nm \times n. 1884년 제임스 조지프 실베스터가 px=xqpx = xq 꼴로 처음 다뤘고, 오늘날 제어이론(관측기 설계·모델 축소), 고유값 문제의 블록 대각화, 영상 복원, 분리가능 편미분방정식의 이산화까지 넓게 퍼져 있다.

핵심 성질은 XX에 대해 완전히 선형이라는 것. 리카티 방정식XSXXSX 같은 이차항이 없어서 원리적으로는 연립 1차방정식일 뿐이다. 문제는 그 “1차방정식”이 미지수 mnmn개짜리라 순진하게 펼치면 감당이 안 된다는 것이고, 실베스터 방정식의 수치해석은 사실상 펼치지 않고 푸는 법에 관한 이야기다.1

2. 유일해 존재 조건[편집]

방정식을 선형사상 S(X)=AX+XB\mathcal{S}(X) = AX + XB로 보면 이야기가 단순해진다. 이 사상은 m×nm\times n 행렬 공간에서 자기 자신으로 가는 선형사상이므로, “모든 CC에 대해 유일해가 존재한다”는 것은 “S\mathcal{S}가 가역이다”, 즉 “S(X)=0\mathcal{S}(X)=0이면 X=0X=0이다”와 같은 말이다.

S\mathcal{S}의 고유값은 크로네커 곱이 바로 알려 준다. vec 항등식으로 펼치면

(InA+BIm)vec(X)=vec(C)\left(I_n \otimes A + B^{\top} \otimes I_m\right)\mathrm{vec}(X) = \mathrm{vec}(C)

이고, 이 계수행렬은 크로네커 합이라 고유값이 정확히 λi(A)+μj(B)\lambda_i(A) + \mu_j(B)다. 따라서 결론은 이렇다.

AX+XB=CAX + XB = C가 모든 CC에 대해 유일해를 갖는 것과 λi+μj0\lambda_i + \mu_j \ne 0(i,j\forall i,j), 즉 σ(A)\sigma(A)σ(B)\sigma(-B)가 만나지 않는 것은 동치다.

AABB 각각이 안정이냐 특이하냐는 아무 상관이 없다. AABB가 둘 다 특이해도 고유값 합만 0을 피하면 유일해가 있고, 반대로 둘 다 아주 멀쩡한 행렬이어도 λi=μj\lambda_i = -\mu_j인 쌍이 하나만 있으면 해가 없거나 무한히 많다. 리아푸노프 방정식에서 조건이 λi+λj0\lambda_i + \lambda_j \ne 0인 것도, 이산시간판(스타인 방정식 XAXB=CX - AXB = C)에서 조건이 λiμj1\lambda_i\mu_j \ne 1인 것도 전부 같은 계산의 변주다.

부호 규약에 한 번은 데이게 되어 있다. 문헌에 따라 AXXB=CAX - XB = C로 쓰기도 하는데, 그러면 조건이 σ(A)σ(B)=\sigma(A)\cap\sigma(B) = \varnothing으로 바뀐다. 부호가 하나 뒤집혔을 뿐인데 “고유값이 겹치면 안 된다”와 “고유값 합이 0이면 안 된다”로 말이 완전히 달라지므로, 논문을 읽을 때 항상 어느 규약인지 먼저 확인해야 한다.

AAB-B의 스펙트럼이 허수축을 사이에 두고 갈라져 있으면 — 예컨대 AABB가 둘 다 후르비츠면 — 해를 적분으로 명시할 수 있다.

X=0eAtCeBtdtX = -\int_{0}^{\infty} e^{At}\, C\, e^{Bt}\, dt

미분해서 대입하면 바로 확인된다. 이 표현이 실무적으로 중요한 이유는, CC가 낮은 계수(C=C1C2C = C_1C_2^{\top}, C1C_1m×rm\times r)이면 피적분함수가 통째로 계수 rr 이하라서 XX도 낮은 계수로 잘 근사된다는 사실을 시사하기 때문이다. 대규모 해법이 전부 이 관찰 위에 서 있다.

3. 왜 크로네커화는 계산에 쓰지 않는가[편집]

위 vec 형태는 유도용으로는 완벽하다. 계산용으로는 재앙이다. m=nm = n일 때 계수행렬이 n2×n2n^2 \times n^2이고, 가우스 소거법으로 풀면

O ⁣((n2)3)=O(n6) 연산,O(n4) 메모리O\!\left((n^2)^3\right) = O(n^6) \ \text{연산}, \qquad O(n^4)\ \text{메모리}

가 든다. n=100n = 100짜리 장난감 문제조차 계수행렬 성분이 10810^8개(배정도로 800 MB)에 연산은 101210^{12}번이다. n=1000n = 1000이면 메모리만 101210^{12}개, 즉 8 TB다.

같은 문제를 바텔스-스튜어트 알고리즘O(n3)O(n^3) 연산에 O(n2)O(n^2) 메모리로 끝낸다. n=1000n = 1000에서 연산량 비가 10910^9배 — 1초와 30년의 차이다. 그래서 실무 코드는 vec 형태를 한 번도 만들지 않는다. 이 행렬은 교과서와 존재 증명 속에서만 산다.

거기에 더해, 크로네커 계는 AA, BB의 희소성을 통째로 잃는다. AA가 삼중대각이어도 IA+BII \otimes A + B^{\top}\otimes IBB의 패턴만큼 블록이 채워지고, 소거를 시작하면 fill-in이 폭발한다. 희소행렬 해법으로 우회하려 해도 얻는 게 별로 없다.

4. 조건수와 분리도 sep[편집]

수치해석적으로 이 문제가 얼마나 어려운지는 분리도(separation)가 잰다.

sep(A,B)=minX0AX+XBFXF=σmin ⁣(InA+BIm)\mathrm{sep}(A, -B) = \min_{X \ne 0} \frac{\|AX + XB\|_F}{\|X\|_F} = \sigma_{\min}\!\left(I_n \otimes A + B^{\top}\otimes I_m\right)

즉 선형사상 S\mathcal{S}의 최소 특이값이다. 섭동 해석을 하면 AA+ΔAA \to A+\Delta A 등에 대해

ΔXFXFAF+BFsep(A,B)(ΔAFAF+)\frac{\|\Delta X\|_F}{\|X\|_F} \lesssim \frac{\|A\|_F + \|B\|_F}{\mathrm{sep}(A,-B)}\left(\frac{\|\Delta A\|_F}{\|A\|_F} + \cdots\right)

꼴이 나오므로, (A+B)/sep(\|A\|+\|B\|)/\mathrm{sep}이 사실상 이 문제의 조건수다.

여기서 반드시 짚어야 할 함정이 있다. sep이 작다는 것과 고유값이 가깝다는 것은 같은 말이 아니다. 고유값이 가까우면 sep이 작은 것은 맞지만, 역은 성립하지 않는다. AA, BB가 강하게 비정규 행렬이면 σ(A)\sigma(A)σ(B)\sigma(-B)가 멀리 떨어져 있어도 sep이 기계 정밀도 근처까지 내려갈 수 있다. 고유값 간격 mini,jλi+μj\min_{i,j}|\lambda_i+\mu_j|는 sep의 상계일 뿐이고, 비정규성이 그 사이를 얼마든지 벌린다. 그래서 “고유값이 잘 분리되어 있으니 안전하다”는 판단은 정규행렬에서나 통한다 — GMRES의 수렴 예측에서 벌어지는 일과 정확히 같은 종류의 배신이다.2

sep을 정확히 계산하려면 mn×mnmn \times mn 행렬의 최소 특이값이 필요하니 실무에서는 하이엄의 노름 추정기 같은 반복법으로 추정한다. LAPACK의 dtrsen이 불변 부분공간의 조건수를 보고할 때 내놓는 값이 바로 이것이다.

5. 어디서 튀어나오는가[편집]

  • 블록 대각화와 불변 부분공간. 블록 상삼각 행렬 [AC0B]\begin{bmatrix} A & C \\ 0 & -B\end{bmatrix}[IX0I]\begin{bmatrix} I & X \\ 0 & I \end{bmatrix}로 상사변환하면 비대각 블록이 AX+XB+CAX + XB + C가 된다. 즉 AX+XB=CAX+XB = -C를 풀면 두 블록이 완전히 분리된다. 슈어 분해에서 대각 블록을 맞바꾸는 재정렬 연산이 국소적으로 이 방정식 하나를 푸는 일이고, 두 블록의 고유값이 가까워 sep이 나빠지면 교환이 불안정해지는 이유도 여기서 나온다.
  • 리아푸노프 방정식. B=AB = A^{\top}, C=QC = -Q인 특수 경우. 대칭해가 나오고 안정성 판정·가제어/가관측 그라미안·H2H_2 노름으로 이어지는 거대한 응용 갈래가 통째로 여기 붙어 있다. 자세한 내용은 해당 문서로.
  • 관측기·제어기 설계. 룬베르거 관측기의 오차 동역학을 원하는 스펙트럼으로 바꾸는 문제나, 출력 궤환으로 고유구조를 배치하는 문제가 실베스터 방정식 형태로 정리된다. 여기서는 XX가 좌표변환 자체다.
  • 모델 축소. 축소차수모델의 유리 보간(rational interpolation)·크릴로프 사영 기법에서, 사영 기저가 만족하는 관계가 실베스터 방정식이다. “전달함수를 어떤 점들에서 정확히 맞춘다”는 조건이 곧 AV+VH=BGAV + VH = -BG 꼴로 쓰인다.
  • 분리가능 PDE. 직사각 격자에서 uxx+uyy=fu_{xx} + u_{yy} = f를 이산화하면 미지량을 행렬 UU로 두었을 때 AU+UA=FA U + U A^{\top} = F가 된다. 격자를 vec으로 펴는 순간 크로네커 곱 라플라시안이 되는 그 구조를, 펴지 않고 그대로 두는 것이다.

6. 사촌들[편집]

같은 “선형 행렬방정식” 가족에 속하는 변형이 여럿이고, 이름이 헷갈리기 쉬우니 한 번에 정리해 둔다.

이름형태유일해 조건
실베스터AX+XB=CAX + XB = Cλi+μj0\lambda_i + \mu_j \ne 0
리아푸노프(연속)AX+XA=QA^{\top}X + XA = -Qλi+λj0\lambda_i + \lambda_j \ne 0
스타인 / 리아푸노프(이산)XAXB=CX - AXB = Cλiμj1\lambda_i \mu_j \ne 1
일반화 실베스터AXB+CXD=EAXB + CXD = E행렬다발 조건
\star-실베스터AX+XB=CAX + X^{\top}B = C별도 이론

일반화판 AXB+CXD=EAXB + CXD = E크로네커 곱으로 펴면 (BA+DC)vec(X)=vec(E)(B^{\top}\otimes A + D^{\top}\otimes C)\mathrm{vec}(X) = \mathrm{vec}(E)가 되고, 유일해 조건이 행렬다발 (A,C)(A,C)(D,B)(D,B)의 일반화 고유값이 겹치지 않는 것으로 바뀐다. 해법도 슈어 분해 대신 QZ 분해(일반화 슈어 분해)를 쓰는 것으로 자연스럽게 확장된다. 서술자 계(descriptor system)의 모델 축소나 특이 행렬다발의 블록 대각화에서 이 형태가 나온다.

XX^{\top}이나 XHX^{H}가 섞이는 \star-실베스터 방정식은 겉보기와 달리 성질이 꽤 다르다. 미지수에 전치가 걸리면 C\mathbb{C} 위에서 선형이 아니라 실선형이 되어, 유일성 조건이 고유값 조건이 아니라 두 행렬다발의 상호 관계로 기술된다. 팰린드로믹 고유값 문제(구조진동·주기계 해석)에서 등장한다.

7. 대규모 희소 문제[편집]

mm, nn10510^5를 넘어가면 슈어 분해 자체가 불가능하고, 애초에 XXmnmn개 성분으로 저장할 수도 없다. 이때는 C=C1C2C = C_1C_2^{\top}가 낮은 계수라는 성질을 이용해 해를 인수 형태 XZ1Z2X \approx Z_1Z_2^{\top}로만 구한다. 근거는 앞서 본 적분 표현과, 두 스펙트럼이 잘 분리되어 있으면 해의 특이값이 지수적으로 감쇠한다는 정리들(Penzl, Beckermann–Townsend 계열)이다.

주요 갈래는 셋이다.

접근매 반복이 하는 일성격
인수형 ADI (ADI 법)이동된 희소계 (A+pkI)v=(A+p_kI)v = \cdots 풀이시프트 파라미터 선택이 성패를 가름
크릴로프 사영크리로프 부분공간법으로 작은 부분공간 생성 후 그 위에서 소형 실베스터 풀이갈레르킨 조건, 잔차로 정지 판정
유리 크릴로프 / 확장 크릴로프A1A^{-1} 또는 (AσI)1(A-\sigma I)^{-1} 곱까지 기저에 포함수렴 훨씬 빠름, 대신 희소 LU 필요

셋 다 안쪽에서 결국 희소 선형계를 반복해 푸는 구조라, 성능이 전처리기 품질에 그대로 종속된다. ADI 시프트를 최적으로 고르는 문제는 두 스펙트럼 영역 위에서 유리함수를 작게/크게 만드는 졸로타레프 문제로 귀착되며, 두 영역을 실구간이나 타원으로 근사해 닫힌 형태 시프트를 쓰는 것이 고전적 처방이다. 요즘은 잔차를 보고 시프트를 그때그때 고르는 적응형 방식이 더 흔하다.3

8. 관련 문서[편집]

9. Footnotes[편집]

  1. 이름이 같아서 헷갈리기 좋은 것이 실베스터 공식(행렬함수를 고유값에서의 값과 프로베니우스 공변량으로 쓰는 공식)과 실베스터 관성 법칙이다. 셋 다 같은 사람 이름이지만 서로 다른 대상이다. 실베스터는 “행렬(matrix)“이라는 용어 자체를 만든 사람이라 이 바닥에 지분이 넘칠 수밖에 없다.

  2. 그래서 논문에서 “sep이 작다”는 말을 만나면 반사적으로 고유값 간격을 떠올리면 안 된다. 두 스펙트럼이 넉넉히 떨어져 있는데도 sep이 기계 정밀도 근처인 예제를 만드는 것은 어렵지 않다. 비정규 행렬은 언제나 우리 뒤통수를 노리고 있다.

  3. 시프트 고르기가 성능의 전부라서, 같은 ADI 코드가 시프트 세트에 따라 20회에 끝나기도 하고 2000회를 돌기도 한다. “일단 돌려” 놓고 잔차가 안 내려간다며 반복 상한만 올리는 것은 여기서 특히 헛수고다. 범인은 항상 시프트다.