[GE] 2. Numerical Solution of the Inverse Kinematics

Jaeyoung Seon·2022년 6월 11일
0

[2022-1H] 게임공학

목록 보기
2/8
post-thumbnail

1. Introduction

nn개의 joint가 있고, θ1,,θn\theta_1,\cdots,\theta_n에서 각 θj\theta_j 값을 joint angle이라고 가정함. (j를 이용하여 θ\theta 값을 순회함)

각 관절의 끝 부분을 end effector라고 함.

IK를 푼다는 것은 end effector가 최대한 target position과 가까워지도록 해를 찾아 자세를 만드는 것임.

몇 가지 정의 사항

  • kk : end effector의 개수. s1,sks_1,\cdots s_k : 각 end-effector의 위치. 보통 root의 xyz 좌표계를 기준으로 정의함.
  • end effector의 위치는 column vector (s1,s2,sk)T(s_1,s_2,\cdots s_k)^T로 정의하고, s\vec{s}라고 씀.
  • target position 역시 column vector (t1,t2,tk)T(t_1,t_2,\cdots t_k)^T로 정의함. tit_ii번째 end effector가 닿아야 하는 target position
  • joint angle 역시 column vector θ=(θ1,,θn)T\theta=(\theta_1,\cdots,\theta_n)^T로 정의함.
  • i번째 end effector의 desired change in position (얼만큼 변해야 하는지)는 ei=tisie_i=t_i-s_i임. (eie_i : end effector의 변화, tit_i : target position, sis_i : current position. goal - current)
  • FK에서는 end effector position을 함수로 나타냄. s=f(θ)\vec{s}=f(\theta) → 각도 θ\theta가 주어졌을 때 end effector position이 어떻게 변해야 한다!
  • 반대로 IK에서는 goal이 주어졌을 때 θ\theta가 어떻게 되어야 하는지를 찾음. sd:θ=f1(sd)\vec{s}_d:\theta=f^{-1}(\vec{s}_d)

기억해야 할 것

  • θ\theta가 주어졌을 때 position이 얻어지는 건 FK, position이 주어졌을 때 θ\theta가 구해지는 것은 IK

(FK : 어떻게 회전해야 저기로 손가락이 가지? vs IK : 손가락이 여기로 가려면 어떻게 회전해야 하지?)

IK를 하기 위해서는 “대략적으로 좋은 것”을 해야 함.

목표에 도달하는 방법의 가지수는 무한하기 때문에 해가 굉장히 많음. 따라서 많은 해 중에 최적의 해를 찾는 good approximation이 중요함.

  • f()f() 함수는 non-linear operator임. → 해가 하나만 존재하지 않음.
  • IK에는 “unreachable”, 즉 어떻게 해도 목표에 닿을 수 없는, 해가 존재하지 않는 경우가 있으며, 이 경우 최대한 목표를 향해 움직이려는 방향으로 해를 가정해야 함.
  • 또한 가장 좋은 solution이 하나가 아닌 경우도 존재함.

⇒ IK가 non-linear operator인 이유

2. Kinetic Angular Variables [3]

기초 of 기초. 2d에서 시작해서 3d로 확장할 것임.

원 위에서 particle이 움직이고 있다고 가정함. 이때 particle의 좌표는 데카르트 좌표 (x,y)(x,y)로 표현할 수 있음.

하지만 x,yx,y만 주어졌다고 해서 점이 “얼마만큼의 속도로 얼만큼” 움직였는지는 알 수 없음. (직선을 따라 움직이는 점의 속도는 쉽게 구할 수 있지만, 원을 따라 움직이는 점의 속도는 쉽게 구할 수 없음)

이 때문에 원 위에서 점의 움직임을 표현하는 간단한 방법은 angle θ\theta를 사용하여 단위 시간 tt에 이 점이 몇 도만큼 회전했는지를 표현하는 것이 가장 직관적임.

여기서는 x축을 기준으로 회전한다고 가정하고, 반지름 RR, 각 θ\theta가 주어졌을 때 x,yx,y는 간단하게 x=Rcosθ,y=Rsinθx=R\cos\theta,y=R\sin\theta로 표현할 수 있음. (ccw를 사용하기 때문에 반시계 방향으로 회전하는 것이 양수, 시계 방향은 음수)


단위 시간 tt에서 점이 움직이는 방향을 각속도 (angular velocity)라고 하며, ω\omega로 표현함.

단위 시간 Δt\Delta t 동안 움직인 각도를 Δθ=θ(t+Δt)θ(t)\Delta\theta=\theta(t+\Delta t)-\theta(t)로 정의함.

이때의 각속도 ω=limΔt0ΔθΔt=dθdt\omega=\lim_{\Delta t\rightarrow0}\frac{\Delta\theta}{\Delta t}=\frac{d\theta}{dt} (속도=거리/시간)

또한 θ(t)=θi+ω(tti)\theta(t)=\theta_i+\omega(t-t_i) 또는 Δθ=ωΔt\Delta\theta=\omega\Delta t (θi\theta_i는 초기 각도, ttit-t_i는 시간 변화 → θi\theta_i에서 출발하여 각속도 ω\omegattit-t_i만큼 이동한 결과가 θ(t)\theta(t)인 것. 시간 x 속도 == 거리 (각))

각속도 ω\omega는 시간에 따라 변화하므로 각가속도 α=limΔt0ΔωΔt=dωdt\alpha=\lim_{\Delta t\rightarrow 0}\frac{\Delta\omega}{\Delta t}=\frac{d\omega}{dt}

따라서 위 식과 비슷하게 ω(t)=ωi+α(tti)\omega(t)=\omega_i+\alpha(t-t_i) 또는 Δω=αΔt\Delta\omega=\alpha\Delta t이고, 이를 θ\theta에 대한 식으로 다시 쓰면

θ(t)=θi+ω(tti)+12α(tti)2\theta(t)=\theta_i+\omega(t-t_i)+\frac{1}{2}\alpha(t-t_i)^2 또는 Δθ=ωiΔt+12α(Δt)2\Delta\theta=\omega_i\Delta t+\frac{1}{2}\alpha(\Delta t)^2

θ\theta를 라디안으로 표현한다면 Δθ\Delta\theta만큼 움직인 호의 길이 s=RΔθs=R|\Delta\theta|임.

각의 변위가 매우 작은 경우 (Δθ0\Delta\theta\rightarrow 0) sΔrs\cong|\Delta\vec{r}|로 표현할 수 있고 (호의 길이가 직선의 길이와 동일해짐), s=RΔθs=R|\Delta\theta|Δt0\Delta t\rightarrow 0을 이용하면 결과적으로 직선 속도 v\vec{v}는 다음과 같이 표현할 수 있음.

v=Rω|\vec{v}|=R|\omega|


Angular velocity in two dimensions [4]

점 p가 회전하고 있는데 VV라는 힘이 작용하고 있는 상황. 이 경우 VVVV_\perpVV_\parallel로 나눠지고, 점 p가 회전하는 힘은 VV_\perp에 의해 결정됨. (VVVV_\perp로 사영 projection한 만큼 힘이 분산됨)

(VV_{\parallel}은 원의 반지름 벡터와 평행한 힘이기 때문에 작용하지 않음)

위 그림에서 rr : 원점 OO와 점 PP를 잇는 위치 벡터, 극좌표 (r,ϕ)(r,\phi). 점 p에는 선속도 v=v+v\text{v}=\text{v}_\parallel+\text{v}_\perp가 작용함. (반지름 벡터와 평행한 v\text{v}_\parallel, 수직인 v\text{v}_\perp)

이때 각속도 ω\omega는 다음과 같음. (v\text{v}_\parallel은 각이 변하지 않으므로 고려하지 않고, v\text{v}_\perp만 고려함)

ω=dϕdt=vr=vsin(θ)r\omega=\frac{d\phi}{dt}=\frac{\text{v}_\perp}{r}=\frac{\text{v}\sin(\theta)}{r} (v=rω\vec{\text{v}}_\perp=r\omega에서 넘어옴)

각속도는 게임의 여러 군데에서 사용됨. (탄막게임, 회전을 이용하는 게임 등)


Practice

Q. 길이 9r의 투석기가 있고, ω\omega의 각속도로 수박을 던짐. 또 길이 r의 투석기가 있는데, 동일한 각속도 ω\omega를 가짐.

이때, 길이 9r 투석기의 선속도 vlargev_{\text{large}}와 길이 r 투석기의 선속도 vsmallv_{\text{small}}은 몇 배 차이가 나는가?

A. vlarge=Rω=9rωv_{\text{large}}=R\omega=9r\omega, vsmall=Rω=rωv_{\text{small}}=R\omega=r\omega. 따라서 vlarge=9vsmallv_{\text{large}} =9v_{\text{small}}


Angular velocity in three dimensions

3차원에서는 회전각을 axis-angle로 표현함.

기존 euler rotation에는 큰 단점이 존재함. → x, y, z축 기준 회전만 가능하다는 것. 꼭 세 가지 축으로만 회전할 수 있는 것이 아닌데 x, y, z축으로 제한시켜버리면 직관적이지 않음. + 짐벌락 문제 등

→ 이러한 문제를 해결하고 임의의 벡터를 기준으로 회전할 수 있게 하기 위해 axis-angle을 사용함. (임의의 축 기준)

axis-angle은 두 가지를 기준으로 표현함.

  • 어떤 벡터를 기준으로 회전할 것인지에 대한 기준 벡터 u\vec{u}
  • 얼만큼 회전할 것인지에 대한 각도 θ\theta

이를 기존 euler rotation으로 표현하려면 우선 x, y, z축으로 내려서 기준 축을 맞춘 다음 rotation, 원래대로 역변환하는 과정을 거쳐야 함. → 너무 복잡하고 비효율적

3차원 공간에서의 각속도 ω\omega (orbital angular velocity)는 다음과 같음.

ω=ω2du=dθdtu=vsin(θ)ru\omega=\omega_{\text{2d}}u=\frac{d\theta}{dt}u=\frac{v\sin(\theta)}{r}u (벡터 형태)

2d에서는 각속도를 quantity (10, 20 등과 같은 값)로 사용한 반면, 3d에서는 벡터로 사용하여 uu 벡터에 곱해줌.

θ\thetarrvv 사이의 각도이기 때문에 각속도 ω\omega는 다음과 같이 다시 쓸 수 있음.

ω=rvsin(θ)urr=r×vr2\omega=\frac{rv\sin(\theta)u}{rr}=\frac{r\times v}{r^2}

벡터 aa, bb가 있을 때, a×ba\times b크기는 두 벡터로 이루어진 평행사변형의 넓이 (quantity)임. 근데 여기에 벡터 n\vec{n}을 곱해주면 외적 벡터가 됨. 이 아이디어를 적용한 것이 위의 식. (rr을 아래 위로 곱해준 이유)

즉, a×b=absin(θ)na\times b=|a||b|\sin(\theta)n이므로 이를 적용하여 다시 쓰면

ω×r=vsin(θ)rursin(90°)v^=v\omega\times r=|\frac{v\sin(\theta)}{r}u||r|\sin(90\degree)\hat{v}_\perp=v_\perp

중요한 공식 (많이 쓰게 됨)

(1) ω=r×vr2\omega=\frac{r\times v}{r^2}

(2) ω×r=v\omega\times r=v_\perp (2d에서는 v=Rω|\vec{v}|=R|\omega|)

3. Angular Momentum

Angular momentum and moment of inertia (각운동량 & 관성 모멘트) [5]

advanced한 개념. 물리적으로 올바르고 정확한 애니메이션을 구현할 때 사용되는 개념.

보통은 관절체를 다룰 때 특정 평면 위에서 다루는데, 이때 각운동량은 벡터의 형태로 다루어야 함.

하지만 IK와 같은 분야에서는 이를 스칼라의 형태로 다룰 수 있음. → why? 어떤 평면을 기준으로 회전할 것인지를 알고 있기 때문.

각운동량은 선운동량의 회전으로 생각할 수 있음. 즉, 선운동량 p\text{p}는 질량 m\text{m}과 선속도 vv에 의해 p=mvp=mv의 형태로 정해짐.

이때 각운동량 LL은 관성 II와 각속도 ω\omega에 의해 L=IωL=I\omega임. (회전체에서는 관성을 고려함)

다양한 물체의 관성 값 → 가장 많이 쓰이는건 long uniform rod→through end

관성 모멘트는 속도를 갖고 움직이는 물체를 변화시킬 때 관성에 대한 저항을 II로 표현한 것. II가 클수록 변화하는데 작용하는 힘이 커짐.

관성 모멘트는 “질량이 회전축에 대해 어느 정도로 분포되어 있는지”에 따라 달라지는데, 보통은 모든 질량이 균일하게 분포하고 있다고 가정함.

4. Jacobian Inverse Methods

4-1. Problem defined by Jacobian

Jacobian JJ는 end effector ss에 대한 모든 체인 시스템에 편미분을 수행한 행렬임.

linear transform은 수식화할 수 있고 비교적 직관적이고 간단함. 근데 non-linear transform이라면?

하나의 수식으로 표현하기 힘들 뿐더러 규칙도 발견되지 않아 어려움. → 이러한 non-linear transform에 사용되는 것이 jacobian

jacobian을 사용하면 도형을 잘게 쪼갤 수 있고, 잘게 쪼갠 곳에서 보면 마치 linear function인 것처럼 보이게 됨.

non-linear function을 linear function으로 근사하는 방식

위 그림에서 실제 모션은 곡선 형태의 움직임. but 이를 국소적인 형태로 보면 (Δs\Delta s가 매우 작음) 모션을 곡선 형태로 정의하나 직선 형태로 정의하나 큰 차이가 없음. → linear function으로 계산 가능


Jacobian JJθ\theta 값에 대한 함수임.

J(θ)ij=(siθj)ijJ(\theta)_{ij}=(\frac{\partial s_i}{\partial \theta_j})_{ij} (i=1,,ki=1,\cdots,k, j=1,,nj=1,\cdots,n) → ii는 end effector를 순회할 때 사용하는 변수, jj는 joint들을 순회할 때 사용하는 변수임.

siθj=vj×(sipj)\frac{\partial s_i}{\partial\theta_j}=v_j\times(s_i-p_j) [1]

sis_i는 end effector의 위치. θj\theta_j는 joint angle이므로 위치 / 각도 == 각속도임. 따라서 위와 같은 식이 나오는 것.

(vjv_j : 회전하고자 하는 축의 unit vector, sis_i : end effector의 위치, pjp_j : joint의 현재 위치)


forward dynamics에서는 s=f(θ)\vec{s}=f(\theta)임. 따라서 이렇게 쓸 수 있음. → dsdtt=J(θ)dθdtt\frac{d\vec{s}}{dt_t}=J(\theta)\frac{d\theta}{dt_t} (이때 J(θ)=dsdθJ(\theta)=\frac{d\vec{s}}{d\theta}, ttt_t : 시간) (양변을 시간에 대해 미분한 결과)

FK에서는 θ,s,t\theta,\vec{s},\vec{t}를 모두 알고 있기 때문에 jacobian J=J(θ)J=J(\theta)를 계산할 수 있음.

즉 end effector의 위치 변화는 Δθ\Delta\theta만큼 업데이트했을 때 ΔsJΔθ\Delta\vec{s}\approx J\Delta\theta임. (위 식과 동일)

IK의 목표는 Δs\Delta\vec{s}가 desired position e=ts\vec{e}=\vec{t}-\vec{s}와 최대한 같아지게 만드는 것.

따라서 FK는 e=JΔθ\vec{e}=J\Delta\theta를 푸는 것이고, IK는 Δθ=J1e\Delta\theta=J^{-1}\vec{e}를 푸는 것임.

4-2. Jacobian Pseudo-inverse

역행렬을 구하는 것이 그렇게 쉽지는 않음. → 행렬의 크기가 커질수록, 또 정사각행렬이 아닌 경우 역행렬 구하기는 훨씬 어려워짐.

→ Jacobian이 정사각행렬이 아니거나 역행렬이 존재하지 않는 경우가 있음.

역행렬이 존재하더라도 singular problem이 발생할 수 있음. (아무리 변화를 시켜도 도달할 수 없는 목표가 있음 → unreachable)

이를 해결하기 위해 Jacobian Pseudo-inverse (J+)J^+)를 정의함. → 반드시 inverse가 나오는 대신, 정확한 inverse가 아니라 근사치를 찾는 것. (따라서 JJ+IJJ^+\fallingdotseq I)

  • 장점 : inverse가 무조건 나온다.
  • 단점 : 정확하지 않다.

⇒ 하지만 정확하지 않더라도 해를 구할 수 있기만 하면 되기 때문에 IK에서는 Jacobian pseudo-inverse를 사용함.

J+J^+는 방정식 e=JΔθ\vec{e}=J\Delta\theta를 “least square” 방식으로 푸는 것임.

least square?

여러 데이터의 집합에서 “오차가 가장 적은” best fit을 찾는 것. 직선의 형태에서 출발하여 2, 3차 함수 등 곡선의 형태로 발전하여 오차를 줄이게 됨.


null space of a matrix

여러 속성 사이의 선형 관계를 정의하기 위해 사용되는 개념. 행렬 A,BA,B에 대해 AB=0AB=0일 때, 행렬 AA의 null space는 항상 zero vector를 포함함. → A0=0A0=0

행렬 (IJ+J)(I-J^+J)JJ의 null space로 projection을 수행한다. → IJ+J=0I-J^+J=0이 되도록 하는 J+J^+를 찾는 것이 pseudo-inverse인 것.

따라서 모든 벡터 ψ\psi에 대해 J(IJ+J)ψ=0J(I-J^+J)\psi=0이 성립함. Δθ=J+e+(IJ+J)ψ\Delta\theta=J^+\vec{e}+(I-J^+J)\psi → IK의 방정식 Δθ=J1Δθ\Delta\theta=J^{-1}\Delta\theta에 더해 (IJ+J)ψ(I-J^+J)\psi가 0이 되는 값을 찾는 것. (0으로 수렴하기 때문에 약간의 에러 값이 추가됨 → 어쨌든 IK를 풀 수 있는 상태가 됨.)

에러 JΔθeJ\Delta\theta-\vec{e}를 최소화하는 방향으로 문제를 품. (약간의 에러는 있을 수 있지만 어쨋든 해를 구할 수 있는 상태를 만드는 것이 목표임 → 이 정도 에러는 눈으로 구별되지 않음.) → optimization problem


Practice

Q. u=(u,v), p=(p,q), x=(x,y)\text{u}=(u, v), \ \text{p}=(p, q), \ \text{x}=(x,y)일 때, T(u)=x\text{T}(\text{u})=\text{x}이고, T(u,v)=(u2v2,2uv)\text{T}(u,v)=(u^2-v^2,2uv)임.

  • u(t)=(t,t2)\text{u}(t)=(t, t^2)t=1t=1에서의 속도는?
  • T(u,v)\text{T}(u,v)의 Jacobian을 구하라.

A.

  • u(t)=(1,2t)\text{u}'(t)=(1,2t)이므로 t=1t=1에서의 속도는 (1,2)(1,2)
  • JT=[T1uT1vT2uT2v]=[2u2v2v2u]J_T=\begin{bmatrix}\frac{\partial T_1}{\partial u} & \frac{\partial T_1}{\partial v} \\ \frac{\partial T_2}{\partial u} & \frac{\partial T_2}{\partial v}\end{bmatrix}=\begin{bmatrix}2u & -2v \\ 2v & 2u\end{bmatrix}

Jacobian matrix의 정의

5. Singular value decomposition

5-1. SVD Introduction

SVD는 팔방미인 → 딥러닝, optimization, pseudo-inverse 등등에 많이 쓰임

SVD는 eigen decomposition & polar decomposition과 관련이 깊음.

행렬의 SVD는 “3개의 행렬로 쪼개는 것”이라고 볼 수 있음. ex) A=BCDA=BCD

그 중 가운데 행렬, 즉 CC에 singular value가 들어가기 때문에 singular value decomposition이라고 부름.

행렬 AA에 벡터 xx를 곱해서 새로운 행렬 AxAx를 만든다고 가정. (변환을 수행할 때 벡터에 행렬을 곱하는 것과 마찬가지)

  • [A]ij[A]_{ij} or aija_{ij}는 행렬 AA에서 “row i와 column j에 있는” 값을 의미함. (a11,a22a_{11}, a_{22} 등)
  • 벡터는 크기와 방향을 가지고 있음.
  • 벡터에 행렬을 곱한다는 것은 “그 벡터에 회전 또는 길이 변화 (stretching)”를 수행한다는 의미임. (affine transform 참고)

matrix factorization (decomposition)

행렬을 쪼개서 원하는 형태로 만드는 것. LU 분해, QR 분해 등 여러 가지 종류의 decomposition이 있음.

5-2. Transpose (전치행렬)

column vector ↔ row vector로 서로 바꾸는 행위. i-th row & j-th column 행렬을 j-th row & i-th column 행렬로 바꾸는 행위

ex) [AT]ij=[A]ji[A^T]_{ij}=[A]_{ji}

ex) C=[123456]C=\begin{bmatrix}1&2\\3&4\\5&6\end{bmatrix}CT=[135246]C^T=\begin{bmatrix}1&3&5\\2&4&6\end{bmatrix}

(AT)T=A(A^T)^T=A임.

중요한 특징 → (AB)T=BTAT(AB)^T=B^TA^T (순서 뒤집혀서 나옴)

또다른 중요한 특징

행렬이 실수가 아닌 경우 ([12j3+i]\begin{bmatrix}1-2j&3+i\end{bmatrix}) → transpose를 취하면 conjugate (켤레복소수)가 됨. (AA^*는 복소수 행렬에 대한 transpose)

[12j3+i]T=[12j3+i]=[1+2j3i]\begin{bmatrix}1-2j&3+i\end{bmatrix}^T=\begin{bmatrix}1-2j&3+i\end{bmatrix}^*=\begin{bmatrix}1+2j\\3-i\end{bmatrix}

5-3. Partitioned Matrix

행렬 B=[135246]B=\begin{bmatrix}1&3&5\\2&4&6\end{bmatrix}이 있을 때, 이를 column matrix로 나타내면 B=[u1u2u3]B=\begin{bmatrix}u_1&u_2&u_3\end{bmatrix}가 됨. → u1=[12],u2=[34],u3=[56]u_1=\begin{bmatrix}1\\2\end{bmatrix}, u_2=\begin{bmatrix}3\\4\end{bmatrix}, u_3=\begin{bmatrix}5\\6\end{bmatrix} 또는 u1T=[12],u2T=[34],u3T=[56]u_1^T=\begin{bmatrix}1&2\end{bmatrix}, u_2^T=\begin{bmatrix}3&4\end{bmatrix}, u_3^T=\begin{bmatrix}5&6\end{bmatrix}

따라서 BT=[123456]B^T=\begin{bmatrix}1&2\\3&4\\5&6\end{bmatrix}, B=[u1Tu2Tu3T]B=\begin{bmatrix}u_1^T\\u_2^T\\u_3^T\end{bmatrix}

행렬을 B=u1x+u2y+u3zB=u_1x+u_2y+u_3z와 같이 나눠서 쓸 때 partitioned matrix를 사용함.

5-4. Eigenvalues & Eigenvectors

어떤 벡터의 방향을 바꾸지 않고 크기만 바꾸고 싶은 경우 스칼라 값을 곱하면 됨.

  • 벡터 uu와 스칼라 값 λ\lambda가 있을 때, 벡터 λu\lambda uuu와 방향은 같지만 크기가 다른 벡터임.
  • 이때 행렬 AAeigenvectorAu=λuAu=\lambda u를 만족하는 non-zero vector uu를 말하고, λ\lambdaeigenvalue라고 함.
  • ex) 행렬 B=[1102]B=\begin{bmatrix}-1&1\\0&-2\end{bmatrix}가 있고, eigenvector u1=[10]u_1=\begin{bmatrix}1\\0\end{bmatrix}, u2=[11]u_2=\begin{bmatrix}-1\\1\end{bmatrix}일 때, 각각의 eigenvalue는?

Bu1=λ1u1Bu_1=\lambda_1u_1, Bu2=λ2u2Bu_2=\lambda_2u_2

[1102][10]=[10]=λ1[10]\begin{bmatrix}-1&1\\0&-2\end{bmatrix}\begin{bmatrix}1\\0\end{bmatrix}=\begin{bmatrix}-1\\0\end{bmatrix}=\lambda_1\begin{bmatrix}1\\0\end{bmatrix}λ1=1\lambda_1=-1

[1102][11]=[22]=λ2[11]\begin{bmatrix}-1&1\\0&-2\end{bmatrix}\begin{bmatrix}-1\\1\end{bmatrix}=\begin{bmatrix}2\\-2\end{bmatrix}=\lambda_2\begin{bmatrix}-1\\1\end{bmatrix}λ2=2\lambda_2=-2

하지만 실전에서는 행렬만 주어지고 eigenvector와 eigenvalue를 구하라고 함. (식은 1개, 미지수는 2개 ㄷㄷ) → 구하는 방법이 있긴 한데, 노가다임. (여기서는 다루지 않음. 컴퓨터가 알아서 계산해 주겠지 뭐)

실전에서는 getEigenVector(), getEigenValue()

5-5. Eigen Decomposition of a matrix [6]

행렬 AA를 eigenvalue와 eigenvector로 분해하는 것은 “행렬 대각화”의 일종임.

SVD는 eigen decomposition의 한 종류임. (eigen decomposition > SVD)

AA를 eigenvalue λ1,λ2,,λk\lambda_1,\lambda_2,\cdots,\lambda_k와 eigenvector X1,X2,,XkX_1,X_2,\cdots,X_k로 나눌 수 있다고 가정하면 다음과 같이 쓸 수 있음

[x11x12x1k],[x21x22x2k],,[xk1xk2xkk]\begin{bmatrix}x_{11}\\x_{12}\\\vdots\\x_{1k}\end{bmatrix}, \begin{bmatrix}x_{21}\\x_{22}\\\vdots\\x_{2k}\end{bmatrix}, \cdots, \begin{bmatrix}x_{k1}\\x_{k2}\\\vdots\\x_{kk}\end{bmatrix}


행렬 PP가 다음과 같이 eigenvector로 이루어진 행렬과 eigenvalue로 이루어진 행렬 DD로 분해된다고 정의함.

P=[X1X2Xk]=[x11x21xk1x12x22xk2x1kx2kxkk]P=\begin{bmatrix}X_1&X_2&\cdots&X_k\end{bmatrix}=\begin{bmatrix}x_{11}&x_{21}&\cdots&x_{k1}\\x_{12}&x_{22}&\cdots&x_{k2}\\\vdots&\vdots&\ddots&\vdots\\x_{1k}&x_{2k}&\cdots&x_{kk}\end{bmatrix} (column vector로 나눈 것)

D=[λ1000λ2000λk]D=\begin{bmatrix}\lambda_1&0&\cdots&0\\0&\lambda_2&\cdots&0\\\vdots&\vdots&\ddots&\vdots\\0&0&\cdots&\lambda_k\end{bmatrix}

위에 나온 식 Au=λuAu=\lambda u를 적용하면 다음과 같이 쓸 수 있음.

AP=A[X1X2Xk]AP=A\begin{bmatrix}X_1&X_2&\cdots&X_k\end{bmatrix}

=[AX1AX2AXk]=\begin{bmatrix}AX_1&AX_2&\cdots&AX_k\end{bmatrix}

=[λ1X1λ2X2λkXk]=\begin{bmatrix}\lambda_1X_1&\lambda_2X_2&\cdots&\lambda_kX_k\end{bmatrix}

=PD=PD

⇒ ⭐AP=PDAP=PD


AP=PDAP=PD에서 AA만 남기고 우항으로 넘기면 A=PDP1A=PDP^{-1}이 됨. → eigen decomposition의 결과

여기서 DD가 eigenvalue로 이루어진 대각행렬이기 때문에 “행렬 대각화”라는 말을 쓴 것.

제약 조건 → AA 항상 정사각행렬이어야 함.

특징)

A2=(PDP1)(PDP1)=PD(P1P)DP1=PD2P1A^2=(PDP^{-1})(PDP^{-1})=PD(P^{-1}P)DP^{-1}=PD^2P^{-1} (DD만 변화함) → (일반화) An=PDnP1A^n=PD^nP^{-1}

역행렬도 마찬가지 → A1=(PDP1)1=PD1P1A^{-1}=(PDP^{-1})^{-1}=PD^{-1}P^{-1}

이때 D1=[λ11000λ21000λk1]D^{-1}=\begin{bmatrix}\lambda_1^{-1}&0&\cdots&0\\0&\lambda_2^{-1}&\cdots&0\\\vdots&\vdots&\ddots&\vdots\\0&0&\cdots&\lambda_k^{-1}\end{bmatrix}λ1,λ2,,λk\lambda_1,\lambda_2,\cdots,\lambda_k의 역수

TIP) 직교행렬인 경우 P1=PTP^{-1}=P^T

(실수에서 정의된 경우에만) PTP=IP^TP=I


Practice

Q. A=[3112]\text{A}=\begin{bmatrix}3 & 1 \\ 1 & 2\end{bmatrix}이고, eigenvector u1=[0.85070.5257], u2=[0.52570.8507]\text{u}_1=\begin{bmatrix}0.8507 \\ 0.5257\end{bmatrix}, \ \text{u}_2=\begin{bmatrix}-0.5257 \\ 0.8507\end{bmatrix}임. eigenvalue가 λ1=3.618, λ2=1.382\lambda_1=3.618, \ \lambda_2=1.382일 때, A=PDPT=?\text{A}=\text{P}\text{D}\text{P}^T=?

A. [0.85070.52570.52570.8507][3.618001.382][0.85070.52570.52570.8507]\begin{bmatrix}0.8507 & -0.5257 \\ 0.5257 & 0.8507\end{bmatrix}\begin{bmatrix}3.618 & 0 \\ 0 & 1.382\end{bmatrix}\begin{bmatrix}0.8507 & 0.5257 \\ -0.5257 & 0.8507\end{bmatrix}


5-6. Singular value

행렬 AA에 대해 ATAA^TA무조건 정사각행렬이 됨. (ex. 3x2 행렬의 경우 전치행렬은 2x3 → 곱하면 3x3 행렬이 됨) 또한 AA의 모든 값이 실수라면 eigenvalue 역시 실수임.

이때 ATAA^TA의 eigenvalue의 square root (\sqrt{})를 singular value라고 함. → λ\sqrt{\lambda}

5-7. Rank

행렬 AA의 rank (계수) → 행렬 AA의 column vector로 나타낼 수 있는 벡터 공간의 수

만약 같은 차원을 나타내는 벡터가 있다면 하나로 취급함.

ex)

따라서 행렬 AA선형 독립인 최대 column 개수가 rank가 됨.

행렬에서 같은 차원을 나타내는 최대 column 개수가 rank와 동일한 경우 full rank라고 함. 아닌 경우 rank-deficient

ex)

행렬 [101231330]\begin{bmatrix}1&0&1\\-2&-3&1\\3&3&0\end{bmatrix}의 rank == 2

첫 번째 column과 두 번째 column은 선형 독립, but (첫 번째 column) - (두 번째 column) = [110]\begin{bmatrix}1\\1\\0\end{bmatrix}이므로 세 번째 column을 정의할 수 있음. 따라서 선형 독립 X

선형 독립?

위 그림의 경우 세 벡터가 모두 다른 평면에 존재하여 다른 두 벡터로 한 벡터를 표현할 수 없을 때 선형 독립이라고 하고,

같은 평면에 있는 벡터가 존재하여 두 벡터로 한 벡터를 표현할 수 있을 때 선형 독립이 아님.

왼쪽 → 3차원, rank = 3 / 오른쪽 → 2차원, rank = 2


Practice

Q. A=[11021102]\text{A}=\begin{bmatrix}1 & 1 & 0 & 2 \\ -1 & -1 & 0 & -2\end{bmatrix}의 rank는?

A. column vector [11]\begin{bmatrix}1 \\ -1\end{bmatrix}이 첫 번째, 두 번째 column에서 서로 같으므로 rank 1, [00]\begin{bmatrix}0 \\ 0\end{bmatrix}은 원점이므로 rank가 늘어나지 않음. [22]\begin{bmatrix}2 \\ -2\end{bmatrix}[11]\begin{bmatrix}1 \\ -1\end{bmatrix}를 2배하여 표현할 수 있으므로 같은 차원에 있음. 따라서 rank = 1


5-8. Singular value decomposition (SVD)

m×nm\times n 행렬 AA의 SVD는 A=UVA=U\sum V^*로 표현할 수 있음.

  • UU : m×mm\times m complex unitary matrix (UU=IU^*U=I)
  • \sum : m×nm\times n (eigen decomposition과 달리 정사각행렬이 아니어도 됨) 행렬. 대각선에 음수가 아닌 실수가 있음.
  • VV : n×nn \times n complex unitary matrix (VV=IV^*V=I)

AA가 실수 행렬인 경우 U,VU,V는 실수 공간에서의 직교 행렬임. 따라서 UVTU\sum V^T로 표현 가능

  • 대각선에 있는 값 σi=ii\sigma_i=\sum_{ii}는 행렬 AAsingular value임. (ATAA^TA의 eigenvalue)
  • 0이 아닌 singular value의 수 == AA의 rank
  • UU의 column을 AA의 left-singular vector, VV의 column을 AA의 right-singular vector라고 함.
  • 이들은 각각 직교 기저 u1,,umu_1,\cdots,u_mv1,,vnv_1,\cdots,v_n을 구성함.

5-9. Geometric meaning of SVD

  • eigenvector는 transform을 통해 어떤 벡터, 혹은 좌표계 자체를 “어떻게" 늘릴 (stretching) 것인지에 대한 의미. (방향 등)
  • eigenvalue는 “얼만큼” 늘릴 것인지에 대한 것.
  • eigenvector는 확대/축소를 해도 span (시점과 종점을 이은 공간)을 벗어나지 않음.

(참고)

(이미지를 클릭하면 영상이 재생됩니다)


SVD의 기하학적 의미

A=UVA=U\sum V^*인데, 이는 곧 벡터에 대한 AA라는 변환을 VV^*\sumUU라는 3가지 단계로 나누겠다는 의미.

T:KnKmT:K^n\rightarrow K^m이라는 선형 변환이 있음.

이때 linear map (AA 행렬에 의해 수행됨)은 RnR^n의 원을 RmR^m의 타원으로 바꾸는 것을 의미함.

0이 아닌 singular value는 변화한 타원의 장축/단축을 의미함. (얼마만큼 stretch)

→ 복잡한 변환을 3개의 덜 복잡한 변환으로 나누어 생각하는 것이 SVD

5-10. SVD calculation (A=UVTA=U\sum V^T)

  • VV의 column은 ATAA^TA의 eigenvector와 같음 (row vector)
  • UU의 column은 AATAA^T의 eigenvector와 같음 (column vector)
  • AAT,ATAAA^T,A^TA의 eigenvalue를 각각 square root하여 나열한 것이 \sum

5-11. SVD example (A=UVTA=U\sum V^T)


Reference

[1] Efficient computation of the jacobian for robot manipulators, 1984

[2] https://www.cuemath.com/data/least-squares/

[3] Motion on a Circle (Or Part of a Circle). (2020, November 6). University of Arkansas. (https://phys.libretexts.org/@go/page/22249)

[4] https://en.wikipedia.org/wiki/Angular_velocity

[5] https://en.wikipedia.org/wiki/Angular_momentum

[6] https://mathworld.wolfram.com/EigenDecomposition.html

[7] https://www.youtube.com/watch?v=PFDu9oVAE-g

profile
재밌는 인생을 살자!

0개의 댓글