이번 포스팅은 지난 강의 Lecture 23-(1)에 이은 미분방정식과 선형대수의 두 번째 이야기이다. 앞선 강의와 이어지는 내용이니만큼 앞 강의 내용을 잘 이해하고 보는 것을 추천한다.
4. 행렬 지수함수(Matrix exponential)
- Diagonalization for differential equation
지난 강의에서 계속 다루어왔던 아래의 식을 다시 생각해보자.
식 (18)의 행렬 A는 벡터 u와 그의 미분값 u'을 연결시키는(couple)역할을 한다. 그런데 우리는 앞선 강의에서 A의 고유값과 고유벡터를 구하여 지수함수(exponential function)로 이루어진 일반해(general solution)를 도출한바 있다. 일반해는 각각의 상수와 고유값, 그리고 고유벡터들로 이루어진 식이며, 선형조합(Linear combination)의 형태로 분리가 되어 있다. 이는 기존의 행렬 A가 고유값/고유벡터에 의해서 분리(uncouple)가 되는 것이다. 이것이 무슨 말인지 이해하기 위해 식 (18)을 고유벡터의 관점에서 다시 정리해보자.
식 (19)는 (18)의 u를 고유벡터와 가중치(weight)곱의 선형조합에 대한 식으로 정리한 것이다. 이를 행렬에 대한 식으로 정의한 것이 (19.1)이다. 여기서 v는 고유값에 대한 지수함수를, S는 고유벡터들로 이루어진 고유벡터행렬을 의미한다. u=Sv를 풀어서 써보면 식 (19.2)와 같이 일반해로 정리되기 때문에 식에 무리는 없다. 이제 식 (19)의 내용을 (18)에 대입하여 A를 분해하여 정리해보자.
(20.1)은 (18)에서 u대신 Sv를 대입하여 정리한 것이다. 이때 S는 상수이기 때문에 좌변의 미분식에서 S는 앞으로 뺄 수 있다. 양변의 좌측에 S의 역함수를 곱해주면 (20.2)와 같이 정리할 수 있는데, 우변의 S^-1AS가 어디서 많이 본 형태이다. 바로 Lecture 22의 행렬의 대각화(diagonalization)에서 배웠던 형태인데, S^-1AS는 Λ(capital lambda)와 같음을 배웠다. 이를 통해 (20.2)를 다시 (20.3)과 같이 정리할 수 있다.
식 (20.3)을 보면 dv/dt=Λv가 되는데, 식 자체의 의미를 놓고 보면 벡터 변수 v를 시간으로 미분한 것은(시간의 변화에 따른 변수 변화량의 비율) 변수에 어떤 상수값을 곱한 것과 같다이다. 이 의미가 성립하려면 v는 식 (19.1)과 같이 지수함수로 정의되어야 한다. 지수함수는 미분한 값이 바로 자기 자신이 되기 때문이다(지수함수의 지수부의 상수가 미분하면 지수함수 앞에 곱해지고 이것이 람다가 되는 꼴). v는 벡터이기 때문에 여러 단일 변수를 동시에 정의할 수 있으며 아래 식은 이에 대한 내용을 정리한 것이다.
식 (21)은 dv/dt=Λv를 2x2크기로 가정하고 식을 풀어 정리한 것이다. 실제로 위의 식과 같이 v에 지수함수를 대입하여 t로 미분해보면 식이 성립함을 알 수 있다. 여기서 중요한 것은 위의 식을 활용하여 미분방정식을 푸는것이다. 우리가 결국 구하고자하는 것은 u이고, v는 식 (19)와 같이 u=Sv의 꼴로 나타냈을 때 S의 각 고유벡터들에 곱해지는 계수(coefficient)에 해당한다. 그러므로 먼저 v의 시간 t에 대한 함수를 정의해보자. 지난 강의 Lecture 23-(1) 에서 배웠듯이 미분방정식의 해(solution)는 지수함수(exponential function)를 도입하여 정의할 수 있다. 단일 변수 대신 다변수에 대해 정의한 행렬에 지수함수를 도입한다는 것이 다를 뿐 기본적인 개념은 같다.
(22.1)은 식 (21)에서 도출한 식이고, 여기서 양변에 dt를 곱해주고 1/v를 곱해준 뒤 적분(integral)을 해준 것이 (22.2)이다. 적분을 하면 (22.3)이 되고, 좌변의 C1을 우변으로 이항하고 자연로그를 없애주기 위해 양변에 지수함수(exponential function)를 취하면 (22.4)가 된다. 절대값은 +, 혹은 -값이 되지만 임의의 값이 될 수 있으므로 +값을 선택하여 e^(c2-c1)자체를 상수 C로 치환하여 (22.5)를 정의할 수 있다.
여기서 우리는 식을 미분방정식의 해를 구하기 위한 함수의 형태로 정리해야 한다. 즉 특정 시간 t에서 해당 미분방정식의의 값이 무엇인지 구할 수 있도록 v의 t에 대한 함수 v(t)의 꼴로 정의해야 한다. v(t)의 꼴로 정의하는 것은 반드시 초기값이 있어야 한다. v에서 t=0인 경우는 지수함수부분이 1이 되므로 상수 c1, c2만 남게 된다. 이 상수 벡터가 v의 초기값이 되고 (22.5)의 C를 대체하여 생각할 수 있다. C를 v의 초기값으로 대체하여 (22.6)과 같이 v(t)를 정의하였다.
v(t)를 구했으니 이제 u(t)를 구해보자. 어떻게 구할 수 있을까? 바로 식 (19)의 u=Sv로부터 구할 수 있다. 아래의 식을 보자.
u=Sv를 u(t)의 꼴로 정리하면 시간의 변수 t를 포함하고 있는 v도 역시 v(t)의 형태로 변해서 (23.2)와 같이 정리할 수 있다. 그 다음 (23.2)의 v(t)를 (23.1)의 v(t)의 우변으로 식을 치환하면 (23.3)과 같이 정리가 된다. u(t)도 마찬가지로 초기값 u(0)가 포함되어야 하는데, u(0)=Sv(0)이고 이로부터 (23.4)와 같이 v(0)를 정의할 수 있다. 여기서 정리한 v(0)를 (23.3)에 대입하여 정리하면 최종적으로 (23.5)와 같이 정리할 수 있다.
정리하자면 우리가 최초에 풀고자 했던 문제는 식 (18)에서 정의한 것과 같이 다변수(multi-variable) 함수에 대한 미분방정식을 푸는 것이었다. 여기서 A는 u와 u의 미분값 사이를 정의하며, 이들을 연결해주는(couple) 역할을 한다. 특정 시간 t에서의 미분방정식의 해(solution)를 구하기 위해선 u를 시간 t에 대한 함수 u(t)의 꼴로 정의를 해줘야 하는데, 이때 사용한 방법이 이 행렬 A의 고유값과 고유벡터를 이용하여 정리하는 것이다. A의 고유값과 고유벡터를 이용하여 대각화를 하고, 여기에 지수함수를 적용하여 식 (23.5)와 같이 미분방정식의 해를 구할 수 있는 식을 만든 것이다. 결국 A의 고유값과 고유벡터를 이용하여 최초의 식 A를 분해하고 식 (19)와 같이 지수함수를 도입하여 다변수 미분방정식의 해(solution)를 구하기 위한 식을 (23.5)와 같이 도출하였다.
그런데 식 (23.5)의 가운데 부분은 다시 아래와같이 정리할 수 있다.
u(t)의 가운데 대각화된 식을 e^At라고 할 수 있는데, 이와 같이 지수함수를 밑으로 하는 형태의 행렬을 행렬지수함수(matrix exponential)라 한다. 뒤이어 알아볼 내용은 식 (24.1)의 행렬지수함수이다.
행렬지수함수를 공부하기 전에 먼저 알아야 할 내용이 있다. 바로 테일러 급수(Taylor series)이다. 인터넷에 테일러 급수에 대한 수많은 자료가 있지만, 매우 중요한 개념이니만큼 간략히 정리하고 넘어가도록 하자.
- Taylor series
테일러 급수(Taylor series)란 어떤 함수를 특정 위치 x=a에서 근사(approximation)하는 방법이다. 테일러 급수의 정의는 아래의 식과 같다.
식 (25)를 풀어서 설명해보자면 n번 미분이 가능한 어떤 함수 f(x)가 있을 때, 이 함수의 특정 위치에서의 값, 즉 x=a에서의 값은 f(x)와 그의 미분값들로 계산한 항들의 무한대의 합으로 표현할 수 있다는 것이다. 이때 각 항은 미분차수값의 factorial로 나누어주어야 한다. 즉 1차 미분이면 1 factorial(1!), 2차 미분이면 2 factorial(2!)로 나누어줘야 한다. 테일러 급수의 식이 왜 이와 같이 표현되었는지를 이해하려면 멱급수(power series)의 전개부터 이해를 해야하는데, 이 부분은 본 포스팅의 범위를 넘어가기 때문에 다음의 동영상을 보는 것을 추천한다(Power series/Euler's great formula).
위에서 설명한 테일러 급수에서 한 가지 중요한 것은 테일러 급수가 근사하는 것은 함수 전체가 아니라 함수의 특정 위치(x=a)에서의 근사만을 의미한다. 물론 어떤 함수의 경우엔 근사하는 위치 주변에서의 값도 잘 맞아떨어질 순 있겠지만, 기본적으로는 x=a에서의 근사를 하기 위함이다. 이와 같은 설명만으론 이해하기 어려울 수 있으므로 예를 한 번 들어보자.
식 (26)의 함수 f(x)는 3차 함수이며 3번 미분이 가능하기 때문에 3차 미분까지만 고려하면 된다. 그 이상은 0이 되기 때문에 의미가 없다. 위의 함수를 x=0일때와 x=-2일때 각각 근사를 해보자. 참고로 x=0일때의 테일러 급수를 매클로린 급수(Maclaurin series)라 한다.
식 (27)은 (26)의 함수 f(x)를 x=0에서 근사한 것이다. (27.1)은 f(x=0)일때의 값과 3차 미분까지의 값들을 나열하였다. (27.2)는 원래의 함수 f(x)를 테일러 급수로 근사한 결과이다. 각 미분항에 곱해지는 (x-a)^n에서 x=0에서의 근사이기 때문에 a에는 0이 대입된다. 정리한 결과 원래의 함수 f(x)가 그대로 복원되고 0을 대입하면 마찬가지로 2가 나오는 것을 볼 수 있다.
이제 x=-2일 때의 근사(approximation)를 수행해보자.
식 (28)은 x=-2에서 테일러 급수를 통해 f(x)를 근사한 모습이다. 마찬가지로 원래의 함수를 잘 복원하는 것을 볼 수 있다.
그런데 여기서 어떤 사람들은 이렇게 생각할 수 있다. "이렇게 그대로 다시 복원할거면 도대체 테일러 급수를 왜 사용하는거지? 그냥 그대로 사용하면 될텐데..". 사실 위에서 예로든 함수 f(x)는 다항식의 유한번의 사칙연산으로 구할 수 있는 대수함수(algebraic function)이고, 따라서 테일러 급수를 통해 원래의 함수를 그대로 복원할 수 있다. 테일러 급수가 진정 빛을 발하는 경우는 대수적 연산만으로는 정확한 값을 구할 수 없고 극한의 과정을 거쳐야만 정확한 값을 구할 수 있는 함수, 이를 테면 지수 함수(exponential function), sin, cos, log등 초월 함수(transcendental function)를 구할 때이다. 초월 함수의 대표적인 테일러 급수의 예로 지수함수(exponential series)와 기하급수(geometric series)에 대한 식을 살펴보자.
식 (29.1)과 (29.2)는 각각 지수함수와 기하 급수를 나타내며, (29.2)의 기하급수는 x의 크기가 1보다 작을 때에만 좌변과 우변의 식이 성립(수렴)한다. 이 외에도 sin, cos, log등의 다양한 초월 함수들이 존재하며 자세한 사항은 테일러 급수 위키(Taylor series)를 참고하기 바란다. 이번 강의를 위해선 일단 위의 두 함수에 대한 테일러 전개(Taylor expansion)를 알면 된다. 참고로 위의 지수함수와 기하급수 두 개의 테일러 급수는 각각 x=0일 때의 매클로린 급수이다.
그렇다면 테일러 급수가 이러한 초월함수에서 힘을 발휘하는 이유는 무엇일까? 그것은 바로 근사(approximation)하는 능력에 있다. 제아무리 성능 좋은 수퍼컴퓨터라해도 무한번의 연산을 수행할 수는 없을 것이다. 따라서 우리는 해당 함수의 값을 원하는 정도의 정밀도가 나올 때까지만 근사하여 구하면 된다. 이때 이러한 근사의 과정을 테일러 급수를 통해서 수행하게 된다. 실제로 공학용 계산기에서 sin함수의 값은 테일러 급수를 통해 어느 정도로 근사된 값이 계산되어 나오는 것이다. 테일러 급수라는 개념이 없었다면 이러한 초월함수들의 계산 과정에 큰 어려움이 있었을 것이다. 근사(approximation)에 대한 이해를 돕기 위해 (29.2)의 기하 급수의 수렴에 대한 예를 살펴보자.
식 (30.1)은 x=1/2일 때의 기하급수의 계산 결과를 나타낸다. 기하급수의 특성상 |x|<1일 때 더해가는 항, 즉 n을 늘려나갈수록 특정 값에 수렴하게 되는데, (30.2)는 두 번째 항까지, (30.3), (30.4)는 각각 세 번째, 네 번째 항까지 더한 결과를 나타낸다. 항이 늘어날 수록 수렴값 2에 가깝게 되는 것을 볼 수 있다. 사실 무한대로 더해야지만 2에 완전히 수렴하겠지만, 무한대로 더하는 것은 불가능하기 때문에 만족할만한 수준의 정밀도를 설정해놓고 그에 상응하는 항까지 더하여 근사를 하면 된다.
이 외에도 테일러 급수는 어떤 복잡한 모델을 단순화 시켜서 풀거나, 미분방정식의 해를 구할 때, 정적분을 구할 때 등 과학과 공학에 있어서 다양하게 활용되는 굉장히 중요한 개념이므로 잘 이해하도록 하자.
- Matrix exponential
앞서 테일러 급수를 간략히 공부했으니 이제 행렬지수함수(Matrix exponential)를 공부해보도록 하자. 미리 말하지만 행렬지수함수를 공부하는 이유는 식 (24)에서와 같이 미분방정식의 해(solution)에 대한 식을 도출하는 과정에서 행렬지수함수가 나오기 때문에 이를 정의하는 법을 알기 위함이다. 우리가 알아볼 행렬지수함수의 형태는 아래의 식과 같다.
우리가 학창시절에 공부하면서 봤던 지수함수의 형태는 대부분 지수함수 e의 지수부에 상수, 혹은 변수가 있는 형태였다. 그런데 식 (31)에 정의된 행렬지수함수는 지수함수의 지수부에 행렬(matrix)이 있는 경우이다. 여기서 행렬 A는 고정된 값이고 t가 기존의 x와 같은 변수이다. 이러한 행렬지수함수는 어떻게 정의할 수 있을까? 그 힌트는 바로 앞서 공부했던 테일러 급수에 있다. 식 (29.1)의 지수함수의 테일러 전개가 (31)의 행렬지수함수에도 그대로 적용되는 것이다. 다만 단일 변수 대신 행렬의 형태로 전개가 되는 것이다. 아래의 식을 살펴보자.
비교해보기 쉽도록 원래의 식 (29.1)과 함께 정리하였다. (29.1)의 1이 (32)에선 단위 행렬 I가 되고, x는 At가 되며 나머지 항들도 같은 패턴을 보이는 것을 알 수 있다. 이때 (29.1)의 1+x 부분이나 (32)의 I + At등의 부분은 사실은 n^0 !, n^1 !로 각각 나누어 준 것이다. 지수함수와 같은 패턴으로 행렬지수함수(Matrix exponential)이 계산된다는건 알겠는데, 정확히 어떤 형태로 나온다는건지 감이 오지 않을 수 있다. 아래의 행렬을 이용하여 행렬지수함수를 계산해 보도록 하자.
식 (33.1)의 행렬 A와 t=1임을 가정하여 행렬지수함수를 계산하였다. 계산 항은 총 10차 까지 계산하였다. 계산 결과값은 (33.3)과 같은 결과를 보인다. 이 값이 과연 맞는 값인지 MATLAB을 통해 확인해보자. MATLAB에서 행렬지수함수를 계산하기 위해선 expm()함수를 이용하면 된다. 아래 그림은 계산 결과와 소스코드이다.
Fig. 1 Matrix Exponential계산 결과와 소스코드. MATLAB 내장함수결과(ans =), 테일러급수 n=10차 계산결과(M =)
Fig. 1의 좌측은 MATLAB을 이용한 계산결과를, 우측은 소스코드를 나타낸다. 좌측 결과에서 ans = 는 MATLAB의 내장함수인 expm()함수를 이용한 계산결과이고 M = 은 테일러급수를 이용하여 10차까지 근사한 결과이다. 근사값이 내장함수값과 거의 일치함을 보이며, 차수를 더해갈 수록 오차는 줄어들게 되고 보다 정확한 근사를 할 수 있다. 이와 같이 행렬지수함수가 테일러 급수로 표현될 수 있음을 증명하였다.
위에서 행렬지수함수를 테일러 급수를 통해 정의하는 법을 배웠는데, 우리는 앞서 지수함수와 더불어 식 (29.2)의 기하급수(geometric series)에 대해서도 배웠다. 기하급수 역시 식 (29.2)의 기존의 형태에 행렬을 대입하여 정리할 수 있다. 아래 식을 보자.
비교를 위해 역시 기존의 식 (29.2)와 함께 정리했다. (29.2)의 1/1-x를 다시 쓰면 (1-x)^-1와 같이 쓸 수 있는데, 식 (34)와 같이 행렬에 대해 정리했을 때 이는 역함수로 해석할 수 있다. 즉 (34)의 (I-At)의 역함수로 생각할 수 있다. 1/(I-At)가 아니다. (29.2)에서 x의 크기가 1보다 작은 경우에만 식이 성립한다고 했다. 마찬가지로 행렬의 형태로 정리한 (34)에서도 이와 같은 조건이 있는데 바로 At의 모든 고유값이 1보다 작은 경우에만 식이 성립한다는 것이다. 이때 중요한것은 행렬 A에 대한 고유값이 아니라 상수 t가 곱해진 At의 모든 고유값이 1보다 작아야 한다는 것이다. 실제 값을 대입하여 계산 결과를 살펴보자.
위의 계산 결과가 실제로 (I-At)의 역함수를 계산했을 때와 같은지 MATLAB을 통해 확인해보자.
Fig. 2 기하급수(Geometric series)의 Matrix Exponential계산 결과와 소스코드. MATLAB 계산결과(ans =), 기하급수 n=100차 계산결과(M =)
Fig. 2의 왼쪽은 MATLAB 계산결과이고 오른쪽은 소스코드이다. V와 D는 각각 At의 고유벡터와 고유값을, mag1과 mag2는 At의 λ1과 λ2의 크기를 각각 의미한다.
다음으로 ans = 는 (I-At)의 역함수를 계산한 결과이며, M은 t=2.2일 때 기하급수를 테일러 급수를 이용하여 계산한 결과이다. 이때 차수는 n=100차 까지 계산하였다. 결과적으로 역함수를 계산한 결과와 급수를 통해 근사한 결과인 M을 비교해보면 거의 같음을 알 수 있다. 물론 차수를 높일 수록 값은 더욱 같아진다. 이로써 기하급수를 행렬지수함수의 형태로 정리한 식 (34)가 조건이 맞을 경우 식이 성립함을 증명하였다. 위의 MATLAB코드에서 t와 A의 값을 바꾸어가며 테스트 해보면 mag1과 mag2가 1이 넘는 경우엔 식이 성립하지 않음을 볼 수 있으니 테스트해보기 바란다.
급수(series)의 관점으로 봤을 땐 식 (32)가 (34)보다는 훨씬 뛰어난 특성을 보인다. (32)의 경우엔 A나 t가 얼마나 크던 간에 항을 더해갈수록 큰 수로 나누어주기 때문에 나중으로 갈 수록 0에 가까워지고 결국 수렴하게 된다. 반면에 (34)의 경우엔 At의 고유값이 1보다 큰 경우 발산할 수 있다. 기하급수의 행렬지수함수를 알아봤지만, 이번 강의에서 우리가 눈여겨봐야할 급수는 바로 식 (32)의 행렬지수함수이다. 이 행렬지수함수를 이용하여 앞에서 다루었던 식 (24)의 행렬지수함수를 증명해보도록 하자. 그 사이에 많은 내용을 다루었기 때문에 편의를 위해서 식 (24)를 다시 써보도록 하자.
우리가 알고자 하는 것은 식 (24.2)가 어떻게 성립하는지, 즉 $Se^{\Lambda t}S^{-1}$를 어떻게 $e^{At}$로 정의할 수 있는지이다. 식 (32)의 행렬지수함수를 이용하여 이를 정의할 수 있다. 아래의 식을 보자.
식 (36.2)는 (36.1)의 행렬지수함수에서 A를 $S \Lambda S^{-1}$로 치환한 것이다. $A=S \Lambda S^{-1}$는 Lecture 22의 행렬의 대각화 부분에서 배웠던 내용이다. (36.2)의 첫 번째 항은 t의 0제곱과 zero factorial이기 때문에 지수함수 e가 사라져서 고유벡터행렬 S만 남은 것이다. (36.2)에서 고유벡터행렬 $S$와 $S^{-1}$는 상수이기 때문에 양쪽으로 뺄 수 있다. 이렇게 S들을 빼서 정리한것이 (36.3)이고, 괄호안에는 고유값 행렬(Λ)에 대한 식이 남는데, 이는 $e^{\Lambda t}$로 볼 수 있다. 따라서 (36.4)와 같이 정의할 수 있고 식 (24.1)을 증명한 것이다.
한 가지 기억해야 할 것은 식 (36.1)은 언제나 성립하는 반면 식 (36.4)는 항상 성립하는 것은 아니라는 것이다. (36.4)가 성립하기위해선 A가 대각화가 가능해야 한다는 것이다. 다시말하면 행렬 A는 n개의 독립(independent)인 고유벡터(eigenvector)를 가져야 한다. 이는 행렬의 대각화를 배울 때 강조했던 내용이며 (36.4)가 성립하기 위해선 반드시 충족시켜야 할 조건이다. 중요한 내용이므로 다시 한 번 강조하니 꼭 기억하도록 하자.
이즈음에서 다시 한 번 생각해보자. 우리가 이렇게 행렬지수함수를 계산하는 법을 공부하는 이유는 무엇인가? 바로 선형연립미분방정식(simultaneous linear differential equation)의 해(solution)를 구하기 위함이며, 그 해는 식 (24)와 같이 u(t)=의 형태이다. u(t)의 식에서 결국 구해야 할 것은 주어진 행렬 A를 기반으로 식 (24.1)의 행렬지수함수를 구하는 것이다. 사실 (36.1)과 같이 테일러급수의 형태로 근사하여 해를 구할 수도 있다. 그러나 계산과정에 있어 비용이 많이 드는 것이 사실이다. 따라서 고유값과 고유벡터를 기반으로 한 (36.4)와 같은 형태로 해를 구하는 것이 훨씬 효율적이다(n개의 독립인 고유벡터가 존재한다는 가정하에서).
행렬 A의 고유값과 고유벡터를 구하는 것은 어렵지 않지만, 문제는고유값 행렬(Λ)에 대한 행렬지수함수이다. 결국 람다에 대한 행렬지수함수를 얻기 위해 (36.1)의 과정을 거쳐야 한다면 (36.4)의 형태가 큰 장점은 없어보인다. 하지만 고유값 행렬(Λ)은 대각행렬이고, 이에 대한 행렬지수함수는 그 계산방법이 훨씬 간단하다는 이점이 있다. 아래의 식을 보자.
식 (37)은 고유값 행렬(Λ)의 행렬지수함수를 나타낸 것이다. 보다시피 고유값만 알면 쉽게 계산할 수 있는 형태이다. 과연 위의 형태대로 계산하는 것이 맞는지 MATLAB을 통해 확인해보자.
Fig. 3 고유값 행렬(Λ)의 행렬지수함수 계산결과(Left)와 MATLAB소스코드(Right)
Fig. 3은 고유값 행렬의 행렬지수함수의 계산 검증 결과를 나타낸다. 왼쪽의 계산결과에서 V와 D는 각각 행렬 A의 고유벡터행렬과 고유값 행렬이고, ans는 MATLAB 내장함수인 expm()으로 고유값 행렬의 행렬지수함수를 계산한결과를, 그리고 expD는 식 (37)의 형태로 직접 계산한 결과이다. 비교해보면 정확히 일치하는 것을 볼 수 있다. 이로써 계산이 훨씬 간편해졌다. 물론 코드가 더 추가되긴 했지만 행렬이 커짐에 따라 이 계산방식이 빛을 발할 것이다. (내장함수인 expm()의 계산 알고리즘이 이미 이와 같은 방식을 사용한다거나, 혹은 더 좋은 방법을 사용한다면 이야기는 달리진다)
정리하자면 최초에 식 (18)과 같이 행렬의 형태로 선형연립미분방정식이 있고 초기값(initial condition)이 주어졌을 때 우리는 식 (24)와 같은 형태로 해를 구할 수 있다. 해를 구하기 위해선 행렬 A의 고유값과 고유벡터를 기반으로 식 (24.2)와 같이 해를 정의할 수 있으며, 이때 고유값 행렬의 행렬지수함수는 (37)의 형태로 계산할 수 있다.
마지막으로 중요한 포인트 한 가지를 다시 짚어보자. 바로 지난 강좌에서 일부 다루었던 내용이지만 중요한 부분이므로 다시 언급하도록 하겠다.
일반적으로 어떤 시스템의 해를 구한다는 것은 시간(t)이 지날수록 그 값이 수렴하는 형태일 것이다. 그렇다면 언제 수렴하는 것일까? 그 해답은 식 (37)의 고유값 행렬에 있다. 식 (24.2)에서 u(t)는 고유벡터행렬 S와 고유값 행렬로 정의되어 있는데, 여기서 고유벡터행렬은 상수이기 때문에 고정되어 있다. 움직이는 것은 t가 곱해져서 영향을 받는 고유값 행렬의 행렬지수함수인데, 식 (37)을 보면 대각원소가 각 고유값에 대한 지수함수로 정의되어 있다. 이것이 의미하는 것은 t가 양의 방향으로 나아간다고 가정했을 경우 해가 수렴하기 위해선 행렬 A의 모든 고유값들이 0보다 작아야 한다는 것이다. 보다 더 정확히 얘기하자면 행렬 A의 모든 실수부(real part)의 고유값이 0보다 작아야지만 시간 t가 지날 수록 해가 안정상태(stable)가 된다. 이는 어떤 시스템 행렬 A를 설계하는데에 있어서 중요한 포인트이므로 잘 기억하도록 하자.
- Test
지금까지 배운 행렬지수함수(Matrix exponential)를 이용한 미분방정식의 해를 구하는 방법이 과연 맞는지 확인해보도록 하자. 지난 강의 Lecture 23-(1) 에서 직접 구했던 해를 이번 강좌에서 배운 행렬지수함수를 이용하여 구해서 결과가 같은지 확인해보자. 아래는 MATLAB을 이용해 두 방식을 테스트 한 결과이다.
Fig. 4 직접 구한 해와 행렬지수함수를 이용한 해의 비교
Fig. 4의 왼쪽이 계산 결과, 오른쪽이 소스코드이다. 왼쪽의 결과에서 A는 행렬, t는 임의의 시간 4.75를 가정하였고 초기값인 v0는 동일하게 [1 0]T으로 설정하였다. Sol=로 표시된 값이 행렬지수함수를 이용하여 계산한 값이고 u=로 표시된 값이 기존에 직접 계산했던 값의 결과이다. 계산 결과가 동일한 것을 볼 수 있다.
5. 마치며(continue)
이번 강좌에선 지난 시간에 이어 미분방정식과 선형대수에 대한 내용을 다루었다. 어떤 시스템의 미분방정식을 세우고 그 해를 u(t)의 꼴로 구하는 과정을 배웠다. 핵심은 결국 시스템 행렬 A의 고유값과 고유벡터에 있으며, 이들을 기반으로한 행렬지수함수를 통해 그 해를 정의할 수 있다. 또한 행렬지수함수를 테일러 급수를 이용하여 정의하였으며, 행렬의 대각화(Diagonalization)를 통해 행렬지수함수를 정의하는 과정을 배웠다. 다음 포스팅에서 미분방정식과 선형대수의 마지막 강좌를 하도록 하겠습니다.