크로네커 곱

1. 개요[편집]

크로네커 곱 Kronecker product
기호	$A \otimes B$
크기	$(m \times n) \otimes (p \times q) = mp \times nq$
혼합곱	$(A \otimes B)(C \otimes D) = AC \otimes BD$
vec 항등식	$\mathrm{vec}(AXB) = (B^{\top} \otimes A)\,\mathrm{vec}(X)$
고유값	$\lambda_i \mu_j$ (전부 곱한 조합)

행렬 두 개를 곱하는 방법은 하나가 아니다. 이건 “각 성분 자리에 상대 행렬을 통째로 박아 넣는” 쪽이다.

크로네커 곱(Kronecker product)은 $m \times n$ 행렬 $A$와 $p \times q$ 행렬 $B$에 대해, $A$의 각 성분 $a_{ij}$를 블록 $a_{ij}B$로 부풀려 만든 $mp \times nq$ 블록행렬이다.

$$A \otimes B = \begin{bmatrix} a_{11}B & a_{12}B & \cdots & a_{1n}B \\ a_{21}B & a_{22}B & \cdots & a_{2n}B \\ \vdots & & \ddots & \vdots \\ a_{m1}B & a_{m2}B & \cdots & a_{mn}B \end{bmatrix}$$

크기 제약이 전혀 없다는 점에 주목할 것. 보통의 행렬곱은 안쪽 차원이 맞아야 하지만 크로네커 곱은 아무 두 행렬에나 정의된다. 이름은 레오폴트 크로네커의 것이지만 정작 그가 처음 쓴 것은 아니고 요한 게오르크 체하푸스(1858)까지 거슬러 올라간다는 것이 정설이다.¹

정의가 왜 하필 이 모양인지는 인덱스를 두 개 쓰다가 하나로 접었다고 보면 곧바로 이해된다. $A$가 인덱스 $i$에, $B$가 인덱스 $k$에 각각 작용하는 연산자라 하자. 두 인덱스 쌍 $(i,k)$를 사전식으로 한 줄 번호 $(i-1)p + k$에 대응시키면, “$i$에는 $A$를, $k$에는 $B$를 동시에 작용시키는” 연산의 행렬 표현이 정확히 $A \otimes B$다. 그래서 성분식이

$$(A \otimes B)_{(i-1)p+k,\ (j-1)q+l} = a_{ij}\,b_{kl}$$

로 깔끔하게 떨어지고, 앞으로 볼 모든 성질이 “각 인덱스에서 따로 일어나는 일”의 번역에 지나지 않게 된다. 인덱스를 접는 순서를 바꾸면 $B \otimes A$가 나온다.

이 연산이 수치해석에서 중요한 이유는 하나로 압축된다. 텐서곱 구조를 가진 문제 — 직사각 격자, 다차원 분리가능 연산자, 복합 양자계 — 는 전부 크로네커 곱으로 쓰이고, 그렇게 쓰인 순간 계산량이 통째로 무너져 내린다. 아래에서 볼 혼합곱 성질과 vec 항등식이 그 무너뜨리는 도구다.

2. 혼합곱 성질 — 거의 모든 것의 원천[편집]

크기가 맞기만 하면

$$(A \otimes B)(C \otimes D) = AC \otimes BD$$

가 성립한다. 혼합곱 성질(mixed product property)이라 부르는 이 한 줄에서 나머지가 거의 다 따라 나온다.

• $(A \otimes B)^{-1} = A^{-1} \otimes B^{-1}$ (둘 다 정칙일 때)
• $(A \otimes B)^{\top} = A^{\top} \otimes B^{\top}$, 켤레전치도 마찬가지
• $A$, $B$가 직교(유니터리)면 $A \otimes B$도 직교(유니터리)
• $A = U_A \Sigma_A V_A^{\top}$, $B = U_B \Sigma_B V_B^{\top}$이면 $A \otimes B = (U_A \otimes U_B)(\Sigma_A \otimes \Sigma_B)(V_A \otimes V_B)^{\top}$ — 즉 특이값 분해가 그대로 크로네커화된다

고유값도 같은 논리로 나온다. 이제 $A$는 $n \times n$, $B$는 $p \times p$ 정사각이라 하자. $Av = \lambda v$, $Bw = \mu w$이면 혼합곱으로

$$(A \otimes B)(v \otimes w) = (Av) \otimes (Bw) = \lambda\mu\, (v \otimes w)$$

이므로 $A \otimes B$의 고유값은 $\lambda_i \mu_j$ 조합 전부($np$개)이고 고유벡터는 $v_i \otimes w_j$다. 같은 이유로 특이값은 $\sigma_i \tau_j$, 계수는 $\mathrm{rank}(A)\,\mathrm{rank}(B)$, 대각합은 $\mathrm{tr}(A)\,\mathrm{tr}(B)$, 행렬식은 $\det(A \otimes B) = (\det A)^{p}(\det B)^{n}$이다. 개별 성분을 하나도 건드리지 않고 $np \times np$ 행렬의 스펙트럼을 전부 알아냈다.

주의할 점 하나. $A \otimes B \ne B \otimes A$가 일반적이다. 다만 둘은 완전 셔플(perfect shuffle) 순열행렬 $P$, $Q$에 대해 $B \otimes A = P(A \otimes B)Q$로 치환 동치라, 고유값·특이값 같은 불변량은 같다. 다시 말해 순서를 바꾸면 “인덱스를 늘어놓는 순서”만 바뀐다.

3. vec 항등식 — 행렬방정식으로 가는 다리[편집]

$\mathrm{vec}(X)$를 $X$의 열을 위에서 아래로 이어 붙인 긴 벡터라 하자(열 우선). 그러면

$$\mathrm{vec}(AXB) = \left(B^{\top} \otimes A\right)\mathrm{vec}(X)$$

가 성립한다. 이 항등식이 하는 일은 명확하다. 미지수가 행렬인 방정식을 미지수가 벡터인 평범한 선형계로 번역해 준다. 크로네커 곱이 수치 선형대수 교과서에서 한 절을 차지하는 이유의 절반이 여기 있다.

즉시 나오는 결과가 실베스터 방정식이다. $AX + XB = C$의 좌변에서 $\mathrm{vec}(AX) = \mathrm{vec}(AXI) = (I \otimes A)\mathrm{vec}(X)$, $\mathrm{vec}(XB) = \mathrm{vec}(IXB) = (B^{\top}\otimes I)\mathrm{vec}(X)$이므로

$$\left(I \otimes A + B^{\top}\otimes I\right)\mathrm{vec}(X) = \mathrm{vec}(C)$$

가 된다. 그리고 잠시 뒤 볼 크로네커 합의 고유값 공식에 의해 이 계수행렬의 고유값이 정확히 $\lambda_i(A) + \mu_j(B)$다. 즉 “$A$와 $-B$가 고유값을 공유하지 않을 때 유일해가 존재한다”는 정리가 행렬식 계산 없이 한 줄로 증명된다. $B = A^{\top}$인 특수한 경우가 리아푸노프 방정식이고, 여기서도 유일성 조건 $\lambda_i + \lambda_j \ne 0$이 같은 자리에서 튀어나온다.

이론적 가치는 이렇게 크지만 계산에는 쓰지 않는다. $n \times n$ 미지 행렬이 $n^2$개 미지수가 되고, 이걸 가우스 소거법으로 풀면 $O(n^6)$ 연산에 $O(n^4)$ 메모리다. 실제로는 바텔스-스튜어트 알고리즘이 같은 문제를 $O(n^3)$에 끝낸다. vec 항등식은 유도용 다리이지 도로가 아니다.

4. 절대 만들지 마라 — 구현 원칙[편집]

같은 경고가 방정식이 아니라 단순 곱셈에도 그대로 적용된다. $(B^{\top}\otimes A)x$를 계산해야 한다고 하자. 정직하게 하면 이렇다.

1. $B^{\top}\otimes A$를 만든다 → $n^2 \times n^2$ 배열, 메모리 $O(n^4)$
2. 벡터에 곱한다 → $O(n^4)$ 연산

vec 항등식을 거꾸로 읽으면 이렇다.

1. $x$를 $n\times n$ 행렬 $X$로 접는다 (reshape, 공짜)
2. $A X B$를 계산한다 → 밀집 행렬곱 두 번, $O(n^3)$
3. 다시 편다 (reshape, 공짜)

연산 $n$배, 메모리 $n^2$배 차이이고, 게다가 2단계가 레벨-3 BLAS라 캐시 효율까지 압도적이다. $n=1000$이면 8 TB 배열을 만드느냐 마느냐의 문제다. 그래서 numpy.kron이나 MATLAB kron은 유도를 확인하는 소형 예제와 교육용 코드에서만 등장해야 하고, 성능이 문제 되는 자리에 나타나면 그건 거의 항상 버그다. “크로네커 곱은 쓰되 만들지는 않는다”가 규칙이다.²

일반화도 자연스럽다. $d$차원 텐서곱 구조 $(A_d \otimes \cdots \otimes A_1)x$는 각 축에 대해 차례로 작은 행렬을 곱하는 것과 같고($d$-모드 곱, 흔히 ttm), 비용이 $O(dn^{d+1})$로 떨어진다. 텐서 계산 라이브러리들이 하는 일의 상당 부분이 이 재배열과 곱셈 순서를 자동으로 고르는 것이다.

5. 크로네커 합과 분리가능 라플라시안[편집]

$A$가 $n \times n$, $B$가 $m \times m$일 때 크로네커 합을

$$A \oplus B = A \otimes I_m + I_n \otimes B$$

로 정의한다. 고유값은 $\lambda_i + \mu_j$, 고유벡터는 여전히 $v_i \otimes w_j$다(혼합곱으로 두 항에 각각 먹여 보면 바로 나온다). 곱이 곱을 만들고 합이 합을 만든다. 덤으로 $e^{A \oplus B} = e^{A} \otimes e^{B}$까지 성립해서 행렬 지수함수 계산이 차원별로 분리된다.

여기서 제일 유명한 응용이 분리가능 2차원 라플라시안이다. 정사각 격자에 유한차분법으로 5점 스텐실을 얹으면, 1차원 2계 미분 행렬 $A = \frac{1}{h^2}\mathrm{tridiag}(1,-2,1)$에 대해 2차원 이산 라플라시안이 정확히

$$L_{2D} = I \otimes A + A \otimes I$$

가 된다(격자점을 열 우선으로 늘어놓았을 때). 그러면 $L_{2D}$의 고유값은 $\lambda_i + \lambda_j$, 고유벡터는 1차원 고유벡터의 텐서곱이다. 1차원 문제의 고유쌍이 $\lambda_k = -\frac{4}{h^2}\sin^2\!\left(\frac{k\pi h}{2}\right)$처럼 해석적으로 알려져 있으므로, $N = n^2$짜리 포아송 방정식을 이산 사인변환 두 번과 대각 나눗셈 한 번으로 $O(N \log N)$에 정확히 풀 수 있다. 고속 포아송 해법(fast Poisson solver)이라 부르는 이것이 고속 푸리에 변환 기반 스펙트럴 방법이나 다중격자 전처리의 정확한 국소 해법으로 지금도 쓰인다.

물론 격자가 비정형이거나 계수가 공간에 따라 변하면 이 구조는 즉시 깨진다. 크로네커 구조는 직사각 격자 + 분리가능 연산자라는 사치스러운 조건의 배당금이다. 그래서 복잡 형상 CFD에서는 이 카드를 못 쓰고, 대신 다중격자법이나 크릴로프 해법으로 간다.

6. 텐서곱 기저와 양자 상태공간[편집]

양자역학에서 복합계의 상태공간은 부분계 상태공간의 **텐서곱**이다. 유한 차원에서 이것을 좌표로 쓰면 그대로 크로네커 곱이다.

$$\dim(\mathcal{H}_A \otimes \mathcal{H}_B) = \dim\mathcal{H}_A \times \dim\mathcal{H}_B$$

차원이 더해지는 게 아니라 곱해진다는 이 한 줄이 양자 시뮬레이션의 모든 고통의 근원이다. 2준위계(큐비트) $N$개면 $2^N$차원이라, $N = 50$에서 상태벡터 하나가 배정도 복소수로 $2^{50} \times 16$바이트 $\approx 16$ PB다. 고전 컴퓨터가 일반 양자계를 직접 시뮬레이션하지 못하는 이유가 이 지수 폭발이며, 파인만이 양자컴퓨터를 제안한 동기이기도 하다.

연산자 쪽도 마찬가지다. $k$번 부분계에만 작용하는 게이트 $U$는 전체 공간에서 $I \otimes \cdots \otimes U \otimes \cdots \otimes I$로 쓰이고, 혼합곱 성질 덕분에 서로 다른 부분계의 게이트는 자동으로 교환한다. 그리고 $\psi = a \otimes b$로 쪼개지지 않는 상태가 바로 양자 얽힘 상태다. 얽힘 정도는 상태를 두 부분계 인덱스로 재배열한 행렬의 특이값 분해(슈미트 분해)로 재고, 이 특이값이 빠르게 감소한다는 사실이 행렬곱 상태·텐서 네트워크 계열 근사가 통하는 이유다.³ 같은 발상이 이징 모형의 전달행렬, 밀도범함수이론 계열의 다체 파동함수 근사에도 그대로 얹혀 있다.

7. 그 밖의 얼굴들[편집]

• 2차원 변환. 2차원 이산 코사인 변환이나 DFT는 $F \otimes F$ 꼴이다. 그래서 2차원 변환을 “행 방향 1차원 변환 + 열 방향 1차원 변환”으로 쪼갤 수 있고, $O(N^2)$이 $O(N\log N)$이 된다. JPEG의 8×8 DCT가 이 분리를 그대로 쓴다.
• 순환행렬 구조. 다차원 합성곱은 블록 순환행렬이고, 이는 크로네커 구조를 가진 DFT 행렬로 동시 대각화된다. 이미지 복원·합성곱 신경망 해석에서 반복적으로 등장한다.
• 최근접 크로네커 근사. 주어진 큰 행렬 $M$에 대해 $\|M - A\otimes B\|_F$를 최소화하는 문제는, $M$의 블록들을 행벡터로 재배열하면 랭크-1 근사 문제가 되어 특이값 분해로 정확히 풀린다(Van Loan–Pitsianis). 전처리기 설계에서 “이 행렬을 텐서곱 구조로 얼마나 잘 흉내 낼 수 있나”를 정량화하는 도구다.
• 통계·기계학습. 다변량 정규분포의 공분산을 $\Sigma_1 \otimes \Sigma_2$로 두는 크로네커 구조 모형은 파라미터 수를 $(nm)^2$에서 $n^2 + m^2$으로 줄인다. 신경망 자연경사법의 K-FAC 근사도 피셔 정보 행렬을 층별 크로네커 곱으로 근사하는 것이 골자다.

8. 관련 문서[편집]

• 실베스터 방정식 · 리아푸노프 방정식 · 바텔스-스튜어트 알고리즘
• 특이값 분해 · 고유값 문제 · 행렬 지수함수 · 행렬식
• 유한차분법 · 포아송 방정식 · 라플라스 방정식
• 고속 푸리에 변환 · 스펙트럴 방법 · 다중격자법
• 양자역학 · 슈뢰딩거 방정식 · 이징 모형
• 희소행렬 · 가우스 소거법 · 자연경사법

9. Footnotes[편집]

1. 이런 “이름 잘못 붙은 정리” 사례가 수학에는 차고 넘친다. 스티글러의 명명 법칙 — “어떤 과학적 발견도 최초 발견자의 이름을 따지 않는다” — 이 스스로에게도 적용된다는 농담까지 있다. 참고로 스티글러는 이 법칙의 최초 발견자로 로버트 머튼을 지목했다. ↩

2. 코드 리뷰에서 kron(이 보이면 일단 크기를 확인하는 습관을 들이는 게 좋다. 8×8짜리면 무해하고, 변수 이름이 n_grid면 그 자리에서 멈춰야 한다. 메모리가 $n^4$로 자라는 코드는 개발 환경의 장난감 격자에서는 아무 문제 없이 돌다가, 실전 격자에서 조용히 OOM으로 죽는다. ↩

3. 슈미트 계수가 급감한다는 것은 곧 “얽힘이 약하다”는 뜻이고, 1차원 갭 있는 기저상태가 면적 법칙을 만족한다는 정리가 이를 보장한다. 반대로 얽힘이 최대인 상태를 만나면 텐서 네트워크도 지수 비용으로 되돌아간다. 공짜 점심은 여기서도 없다. ↩