선형 대수 기초 에서

선형 대수 기초 소개

"선형대수"란 무엇일까요?

위키백과나 교과서, 그리고 3Blue1Brown 같은 훌륭한 웹사이트에 기술적인 정의가 많이 있지만, 우리 목적에는 좀 더 편하게 설명하는 게 좋겠어요.

Julia (그리고 Python) 메인테이너들은 이런 걸 아주 좋아해요. Exercism 경영진은 그렇게까지는 아니고요.

이 개념에서는 거대한 주제 중 더 단순한 부분만 다룰 거예요. 하지만 미리 말해두자면, 이건 필연적으로 꽤 수학적인 개념이에요.

벡터란 무엇일까요?

누구에게 묻느냐, 그리고 다소 추상적인 것을 어떻게 시각화하고 싶으냐에 따라 달라요.

• Julia에서 Vector{T}는 그냥 Array{T, 1}의 별칭이고, eltype인 T는 무엇이든 될 수 있어요.
• (많은) 물리학 분야에서 벡터는 N차원 공간에서 length와 direction을 가진 화살표지만, 고정된 위치는 없어요.
• 선형대수에서 물리학자의 화살표는 공간에 고정되어 꼬리가 원점에 있어요. 그러면 화살표의 몸통은 불필요해지므로, 벡터를 N차원 공간의 point로 표현할 수 있고, N개의 원소가 N개의 각 축을 따라 원점에서 떨어진 거리를 나타내요.

지금 목적에서는 문자열이나 문자의 벡터는 무시해요. 이 개념에서는 숫자 타입, 즉 Int, Float, Complex를 다룰 거예요.

행렬이란 무엇일까요?

이번에도 다양한 답이 있을 거예요.

• 숫자를 직사각형 2차원 배열로 늘어놓은 것 (삐죽 튀어나온 부분 없이). 정사각 행렬은 흔한 특수한 경우예요.
• 열 벡터를 나란히 쌓아 놓은 linear combination.
• 다른 입력 타입에 function을 적용하듯, 벡터에 적용할 수 있는 transformation.

정사각 행렬 중에는 특별한 이름을 가질 만큼 흔한 것들이 있어요.

대각 행렬: 0이 아닌 항목이 모두 main diagonal(왼쪽 위에서 오른쪽 아래)에 있어요.

julia> [1 0 0; 0 2 0; 0 0 3]
3×3 Matrix{Int64}:
 1  0  0
 0  2  0
 0  0  3

# a convenient shortcut
julia> diagm(1:3)
3×3 Matrix{Int64}:
 1  0  0
 0  2  0
 0  0  3

단위 행렬: 대각선에 1만 있는 대각 행렬이에요 (이유는 나중에 더 분명해질 거예요). 흔히 I로 줄여 써요.

julia> Matrix{Float64}(I, 3, 3)
3×3 Matrix{Float64}:
 1.0  0.0  0.0
 0.0  1.0  0.0
 0.0  0.0  1.0

상삼각 행렬: 대각선 위와 대각선에 0이 아닌 값이 있고, 그 아래는 0이에요. 하삼각 행렬 은 짐작하겠죠.

행렬의 전치

행과 열을 정말로 맞바꾸고 싶다면, permutedims() 함수가 그렇게 해서 새로운 행렬을 만들어 줘요.

이건 복사를 수반해서, 큰 행렬을 다룰 때는 느리고 메모리를 많이 먹어요.

선형대수 목적에는 transpose() 함수가 더 유용해요. 원래 행렬을 감싸는 지연 래퍼를 빠르게 만들어 주거든요.

더 유용한 건 adjoint() 함수로, 복소수의 허수부 부호까지 뒤집어 줘요 (왜 이게 그렇게 중요한지 궁금하다면, 한 가지 짧은 답은 양자역학이에요). 이 연산은 워낙 흔해서 변수 이름에 작은따옴표 '를 붙이는 것만으로 adjoint를 만들 수 있어요.

julia> m
2×3 Matrix{Int64}:
 1  2  3
 4  5  6

# new, full copy
julia> permutedims(m)
3×2 Matrix{Int64}:
 1  4
 2  5
 3  6

# lazy version
julia> transpose(m)
3×2 transpose(::Matrix{Int64}) with eltype Int64:
 1  4
 2  5
 3  6

julia> mc = [1+2im 2+3im; 3+2im 1+2im]
2×2 Matrix{Complex{Int64}}:
 1+2im  2+3im
 3+2im  1+2im

# lazy conjugate transpose
julia> adjoint(mc)
2×2 adjoint(::Matrix{Complex{Int64}}) with eltype Complex{Int64}:
 1-2im  3-2im
 2-3im  1-2im

# syntactic sugar for adjoint
julia> mc'
2×2 adjoint(::Matrix{Complex{Int64}}) with eltype Complex{Int64}:
 1-2im  3-2im
 2-3im  1-2im

곱셈

이 문서의 나머지 부분에는 LinearAlgebra 모듈의 기능이 나와요. using LinearAlgebra 한 줄이면 네임스페이스로 불러올 수 있지만, 예제마다 반복해서 쓰지는 않을게요 (시각적으로 너무 어수선하니까요).

벡터, 요소별 곱

벡터 연산 개념에서 다뤘어요. 연산자는 .*이고, 입력 벡터끼리 짝을 지어 연산해서 입력과 같은 크기와 타입의 결과를 내요.

julia> [1, 2] .* [3, 4]
2-element Vector{Int64}:
 3
 8

벡터, 내적

이 아주 흔한 연산은 교과서에서 u ⋅ v로 쓰고 내적이라고 불러요.

내적은 요소별 곱의 _합_과 같아요.

두 벡터는 보통의 * 연산자로 곱할 수 있지만, 왼쪽 벡터가 adjoint일 때만 그래요. 편리하게 u' * v라고 쓰죠. 이건 행렬 곱셈을 다루는 뒷부분에서 더 분명해질 거예요.

위에 찍는(가운데) 점은 Julia에서 (\cdot 탭으로 입력) dot() 함수의 문법 설탕으로 쓸 수 있어요. adjoint를 지정할 필요는 없어요. 이 부분은 자동으로 처리되거든요.

julia> using LinearAlgebra

# with element-wise syntax
julia> sum([1, 2] .* [3, 4])
11

# with u' * v syntax
julia> [1, 2]' * [3, 4]
11

# with dot()
julia> dot([1, 2], [3, 4])
11

# with \cdot syntax
julia> [1, 2] ⋅ [3, 4]
11

복소수 벡터의 경우 왼쪽 벡터가 켤레여야 해요 (허수부의 부호가 뒤집힌). dot 함수는 이걸 자동으로 해주고, u' 문법은 명시적으로 하고, sum(u .* v)는 u와 v가 복소수면 실패해요.

벡터, 외적

예전에 전자기학 강의를 들은 분이라면 향수를 느낄 만한 내용이에요! 돌림힘이나 각운동량 벡터를 계산해 본 엔지니어들에게도요.

내적이 두 벡터를 scalar로 바꾸는 반면, 외적은 두 3-벡터를 세 번째 3-벡터로 바꿔요.

벡터 공간의 기하학적 표현에서, 새 벡터는 두 입력 벡터를 포함하는 평면에 _수직_이에요. 두 입력이 (부호를 빼고) 평행하면 평면을 정의하지 못하므로, 결과는 [0, 0, 0]이 돼요.

순서가 중요해요: u × v == -(v × u). 이게 유명한 오른손 법칙이고, 우리 중 많은 사람이 엄지와 두 손가락을 공간에서 이리저리 돌려 보며 멍하니 바라보게 만든 것이죠 (흔히 어리둥절한 표정으로요).

# with \times syntax
julia> [1, 2, 3] × [3, 4, 5]
3-element Vector{Int64}:
 -2
  4
 -2
# with cross()
julia> cross([1, 2, 3], [3, 4, 5])
3-element Vector{Int64}:
 -2
  4
 -2

이 연산은 길이가 3인 벡터로 제한돼요 (직교하는 x, y, z 축을 가진 유클리드 공간과 같아요). 수학자들은 다른 차원에 대해서도 기꺼이 이 연산을 정의하는데, 그 결과는 이해하기 어려운 용어 잔치(고차 텐서, 쐐기곱, 외부곱, 다중 벡터...)예요. 대부분은 그쯤에서 도망가거나 재빨리 화제를 돌리죠!

노름

벡터는 얼마나 "클까요"?

norm은 벡터를 적절한 스칼라로 줄여서 이걸 잡아내려는 시도예요.

온갖 노름이 정의되어 있지만, 단연 가장 흔한 것은 2-노름이고, √(v ⋅ v)예요.

이 제곱평균제곱근 연산은 (N차원 공간에서) 원점으로부터의 _피타고라스 거리_예요. 벡터를 꼬리가 원점에 있는 화살표로 시각화하면, 2-노름은 화살표의 _길이_예요.

어떤 p-노름이든 p를 두 번째 인자로 넘겨서 계산할 수 있어요. 1-노름이 가끔 유용한데, 항목들의 절댓값의 합일 뿐이에요 (그래서 계산이 아주 빠르고 쉬워요).

# defaults to the 2-norm
julia> norm([1, 2, 3])
3.7416573867739413
# the 1-norm
julia> norm([1, -2, 3], 1)
6.0

행렬 곱셈

세상의 모든 응용 수학자가 지난 80년 동안 온갖 계산을 일련의 행렬 곱셈으로 바꾸는 데 많은 시간을 쏟은 것처럼 보일 때가 있어요.

이건 컴퓨터가 아주 잘하는 연산이에요:

• 반복이 아주 많아요.
• 효율적으로 병렬화할 수 있어요.
• 더 빠르게 만들기 위해 많은 특수 하드웨어가 개발됐어요. 1970년대 8천만 달러짜리 Cray-1부터, 아마 노트북에 통합돼 있을 GPU까지요.

세부 내용은 꽤 단순해요. 다만 처음 보면 그리 직관적이지는 않아요.

행렬 A가 벡터 v를 곱해 결과 w를 낸다고 해봐요 (선형대수 관례상 행렬은 대문자, 벡터는 소문자로 써요).

A의 첫 번째 행을 v와 내적하면 w의 첫 번째 원소가 나오고, 두 번째 행은 두 번째 원소를 주는 식으로 계속 내려가요.

julia> A = [1 2; 3 4]
2×2 Matrix{Int64}:
 1  2
 3  4
julia> v = [5, 6]
2-element Vector{Int64}:
 5
 6
julia> A * v
2-element Vector{Int64}:
 17  # equals [1, 2] ⋅ [5, 6]
 39  # equals [3, 4] ⋅ [5, 6]

행렬 * 행렬 곱셈은 이걸 오른쪽 행렬의 열 전체로 확장해요.

C = A * B에서, A가 B의 각 열을 곱해 C의 대응하는 열을 낸다고 볼 수 있어요. 일련의 행렬-벡터 곱셈이죠.

같은 말로, A의 첫 번째 행을 B의 각 열과 내적하면 C의 첫 번째 행이 나오고, 두 번째 행은 두 번째 행을 주는 식으로 계속 내려간다고 할 수 있어요.

이런 머릿속 표현 방식이 여러 가지가 있는데, 관심 있는 학습자는 이걸 다루는 MIT 강의를 통째로 볼 수 있어요.

julia> A
2×2 Matrix{Int64}:
 1  2
 3  4

julia> B = [5 7; 6 8]
2×2 Matrix{Int64}:
 5  7
 6  8

# left column is the same as A*v previously
julia> A * B
2×2 Matrix{Int64}:
 17  23
 39  53

julia> B * A
2×2 Matrix{Int64}:
 26  38
 30  44

위 예에서 보듯, 행렬 곱셈은 교환 법칙이 성립하지 않아요: A*B와 B*A 사이에는 단순한 관계가 없어요.

행렬 곱셈은 처음 접하는 사람이라면 글만 읽어서는 시각화하기 어려울 거예요. YouTube에는 이걸 그래픽으로 보여주는 영상이 많으니, "행렬 곱셈"을 검색해서 원하는 스타일과 세부 수준, 언어의 영상을 골라 보세요.

차원

두 벡터의 내적은 두 길이가 같아야 성립해요.

이를 확장하면, 행렬 곱셈에서는 왼쪽 행렬의 _열 개수_가 오른쪽 행렬의 _행 개수_와 같아야 해요.

size(A)가 내놓는 것처럼 크기를 (nrows, ncols) 튜플로 나타내면 (a, b) * (b, c) -> (a, c)가 돼요. 여기서 "안쪽" 차원인 b와 b는 내적을 취하는 것과 호환돼요. 여기서 "바깥쪽" 차원인 a와 c는 결과의 차원을 결정해요.

직사각 행렬을 쓴 예:

julia> D = reshape(1:6, 2, 3)
2×3 reshape(::UnitRange{Int64}, 2, 3) with eltype Int64:
 1  3  5
 2  4  6

julia> E = reshape(1:12, 3, 4)
3×4 reshape(::UnitRange{Int64}, 3, 4) with eltype Int64:
 1  4  7  10
 2  5  8  11
 3  6  9  12

julia> D * E
2×4 Matrix{Int64}:
 22  49   76  103
 28  64  100  136

# (3, 4) * (2, 3) not possible
julia> E * D
ERROR: DimensionMismatch: matrix A has axes (Base.OneTo(3),Base.OneTo(4)), matrix B has axes (Base.OneTo(2),Base.OneTo(3))

벡터 쌍을 행렬 방식으로 곱하는 데는 두 가지 가능성이 있어요.

관례적으로는 내적과 같은 것으로 u' * v를 써요. 차원은 (1, 3) * (3, 1)이고, Julia는 (1, 1)짜리 결과를 스칼라로 단순화해요 (예를 들어 R과는 달라요).

또는 u * v'를 써서, 차원 (3, 1) * (1, 3) -> (3, 3)으로 가능한 모든 원소 쌍을 곱해 벡터를 행렬로 확장할 수 있어요.

julia> u = [1, 2]
2-element Vector{Int64}:
 1
 2
julia> v = [3, 4]
2-element Vector{Int64}:
 3
 4
julia> u' * v
11
julia> u * v'
2×2 Matrix{Int64}:
 3  4
 6  8

u' * v는 가끔 _내적_이라고도 불려요.

(덜 흔한) u * v'는 그에 대응해 _외적_이에요. 이건 텐서곱과 관련이 있지만, 그건 우리 범위를 훨씬 벗어나요.

회전

특히 흔한 행렬 곱셈 유형 중 하나는 회전 행렬이에요.

2차원에서는 벡터를 반시계 방향으로 θ 라디안만큼 회전시키는 비교적 단순한 행렬이 있어요.

julia> rot2d(θ, vec) = [cos(θ) -sin(θ); sin(θ) cos(θ)] * vec
rot2d (generic function with 1 method)

# unit vector in the x direction
julia> i_hat = [1, 0]
2-element Vector{Int64}:
 1
 0

# rotate 45 degrees
julia> rot2d(π/4, i_hat)
2-element Vector{Float64}:
 0.7071067811865476
 0.7071067811865475

# rotate 90 degrees -> unit vector in the y direction, j_hat
julia> rot2d(π/2, i_hat)
2-element Vector{Float64}:
 6.123233995736766e-17  # zero, within numerical error
 1.0

3차원 회전에는 두 각도를 가진 더 복잡한 행렬이 필요해요. 이건 LaTeX를 많이 넣지 않고는 Markdown 문서에 보여주기 어려우니, 공식은 위키백과에서 확인해 봐요.

OpenGL과 그 후속작 같은 3차원 그래픽에서는 동차 좌표로 작업하는 게 관례예요.

각 꼭짓점은 4-벡터(또는 동등하게 행렬의 열)로 표현돼요. (x, y, z)에 있는 점이라면 [x, y, z, 1.0]이에요.

이렇게 하면 더 정교한 변환 행렬을 쓸 수 있어요. 회전은 여전히 A[1:3, 1:3]에 있고, 평행 이동은 A[1:3, 4]의 [Δx, Δy, Δz], 크기 조정은 대각선에 있고, 비틀기나 원근법 등 다른 가능성도 있어요.

성능

전 세계 컴퓨터과학과에는 행렬 곱셈 알고리즘을 조금씩 점진적으로 개선한 내용을 한 장 이상 담은 박사 학위 논문이 많아요. 이건 정말 중요하고, 여러 기관이 자기 이익을 위해 그 연구에 기꺼이 자금을 대요.

성능은 깊이 들어갈 수 없는 큰 주제지만, 예를 보여주기 위해 여러 크기의 무작위 행렬을 곱해 볼 수 있어요.

아래 코드는 대충 만든 어림값이에요 (더 나은 방법은 BenchmarkTools.jl을 써 봐요).

사용한 시스템은 430달러짜리 소형 (미국) PC였어요: Ryzen 9 프로세서, 32GB RAM, Linux Mint 22.1, Julia 1.11.6.

julia> mmul(A) = A * A
mmul (generic function with 1 method)
julia> A = rand(Float64, 1_000, 1_000);
julia> @time mmul(A);
  0.011055 seconds (3 allocations: 7.629 MiB)
julia> A = rand(Float64, 10_000, 10_000);
julia> @time mmul(A);
  7.236284 seconds (1.61 k allocations: 763.027 MiB, 17.92% gc time, 0.15% compilation time)

대략 백만 원소짜리 행렬 한 쌍은 몇 밀리초가 걸렸고, 1억 원소짜리 행렬은 몇 초가 걸렸어요. 몇 자릿수를 더 키우려면 더 좋은 하드웨어가 필요하겠죠...