FDM #5, Heat Equation의 풀이(1)

이번 글은

2022.07.25 - [수학의 재미/아름다운 이론] - 미분방정식을 연립방정식으로, FDM #4 (2계 상미분방정식)

미분방정식을 연립방정식으로, FDM #4 (2계 상미분방정식)

에서 이어집니다.

 

Heat equation(열방정식)은 자연 현상에서 열이 이동하는 모습을 설명하는 가장 유명한 편미분 방정식 중 하나입니다.

이 글의 목적은 Heat equation을 FDM으로 풀어보는 것입니다.

 

Heat Equation 이란?

 

다양한 Heat equation 중 가장 간단하고 의미있는 형태는 다음과 같습니다(위키피디아 발췌)

 

Heat Equation

$u(t,x)$가
$$ \frac{\partial u}{\partial t} = \alpha \frac{\partial^2 u}{\partial x^2} \tag{1}$$
형태의 편미분 방정식을 만족할 때, 이를 heqt equation 이라 합니다. 여기서

○ $t$ : 시점을 나타내는 변수, $t\geq 0$
○ $x$ : 변위를 나타내는 변수(공간속의 좌표) $x \in [0,L]$
○ $\alpha$: 상수

○ 초기조건(initial condition) :  $u(0,x)= f(x)$, $\forall x \in [0,L]$
○ 경계조건(boundary condition) : $u(t,0) = 0,  u(t,L)=0$

위키피디아에 수록된, 열이 막대기안에서 전달되는 과정

즉, 초기시점에 공간속에 열이 $f(x)$형태로 주어져 있을 때, 이 열이 시간이 지나며 어떻게 퍼져나가는지를 말해주는 방정식이죠.(막대기의 끝단들에서는 열의 양이 0인 상태)

일반적으로 이러한 형태의 편미분 방정식은 

• 편미분 방정식
• 초기조건(initial condition): $t=0$일 때의 확정된 조건
• 경계조건(boundary condition) : 공간의 경계에서 확정된 조건

이 세 요인으로 해가 결정됩니다.

 

그림으로 본 Heat Equation의 구조

 

 

FDM 구조로 바꾸면?

식(1)의 열방정식을 해결하기 위해서는 식(1) 자체를 이산화(discretization) 시켜야 합니다. 다음의 방법대로 하면 되죠.

 

열방정식의 이산화 과정

1. 시간 $t$ 축을 유한개의 균등 격자로 나눔
  - 식(1)에는 초기시점 $t=0$은 있지만, 끝시점은 없죠. 따라서 관심있는 충분히 큰 끝시점 $T$를 하나 정하고,
  $$ 0= t_0 < t_1 <t_2 <\cdots < t_{N-1} < t_N =T$$
    로 균등하게 나눔
  - 각 시점의 차이는 $k$로 일정함. 즉, $t_{n+1}-t_n =k$

2. 변위 변수 $x$도 유한개의 균등 격자로 나눔
 - 식(1)의 $x$가 움직이는 구간을 $[x_m, x_M]$이라 했을 때,
$$ x_m =x_0 < x_1 <x_2 < \cdots < x_{J-1} < x_J = x_M$$
   으로 균등하게 나눔
- 각 변위의 차이는 $h$로 일정함. 즉, $u_{j+1}-u_j = h$ 

그림으로 설명하면 이렇습니다.

 

시간축과 변위축을 균등하게 잘라 격자판을 만듬

이렇게 격자를 설계하는 이유는

 

 

$U(t,x)$를 모두 구하고 싶은데, 그러기 힘든 경우엔 $U$를 대표하는 점들을 설정하여 그 점들에서의 값만 알아보자.

입니다. 따라서 격자를 잘게 잘게 나눌수록 $u$를 더 자세하게 구하는 효과가 있겠죠. 하지만, 반대급부로 격자가 많아질수록, 계산량이 많아져 원하는 결과가 신속하게 나오지 않게 될 수도 있습니다.

 

자, 그럼 예시를 하나 들어보겠습니다.

 

 

Heat Equation 예시

다음과 같은 문제가 주어져 있다고 합시다.

  $$u_t(t,x) = u_{xx}(t,x),  0<x<1, t>0,$$

○ 초기조건: $ u(0,x) = \sin(\pi x) $
○ 경계조건: $ u(t,0)= u(t,1) =0 $
 
※ 참고 이 방정식의 exact solution 은
$$u(t,x)= \sin(\pi x)\exp (-\pi^2 t)\tag{2}$$

이 방정식의 exact solution이 식(2)처럼 나오는 이유를 스케치만 해보면 다음과 같습니다. 우선 $u(t,x)$가 변수분리되는 함수라고 생각해 봅시다. 즉,

$$ u(t,x) = f(t) g(x)$$

이렇게요. 그러면 위 예제는 다음과 같이 바뀝니다.

 

$$ f'(t)g(x) = f(t) g''(x) ~,~ f(0)g(x)=\sin(\pi x)~, ~f(t)g(0)=f(t)g(1)= 0$$

 

$\frac{1}{f(0)}=c$라고 해보죠.  그러면 우선 $g(x) = c\sin(\pi x)$ 입니다. 이것을 첫째 식에 넣으면

$$ cf'(t) \sin(\pi x) = c^2 f(t) (-\pi^2 \sin(\pi x))$$

을 만족하므로 $f'(t) =-\pi^2 f(t)$ 가 됩니다. 따라서

$$ \left(\ln(f(t)) \right)' = -\pi^2 $$, 즉  어떤 상수 $a$가 있어

$$ f(t) = a e^{-\pi^2 t}$$

의 형태이죠. $t=0$을 대입하면 $f(0)=a$가 나옵니다. 이 값이 $1/c$입니다. 따라서

$$ u(t,x) =f(t) g(x) = ac  \sin(\pi x)\exp (-\pi^2 t) $$

이므로 exact solution을 얻을 수 있죠.(이건 정말 방만하게 푼 해법이고, 자세한 풀이는 어려워서 생략합니다.)

 

 

이제 이 문제를 FDM으로 풀어보도록 하겠습니다. 푸는 방법에는 크게 두 가지 방법이 있는데요. 다음 글에서 자세하게 언급하도록 하겠습니다.