[수치해석] Newmark method(뉴마크 기법)

이번 포스팅에서는 수치 적분 기법 중 하나인

 

Newmark method(뉴마크 기법)에 대해 알아보도록 하겠습니다.

 

 Newmark method(뉴마크 기법) 개요

Newmark method는 대표적인 implicit method 중 하나인데요.

 

Nathan Mortimore Newmark(1910 - 1981)라는 사람에 의해

 

1959년 발표된 논문에 소개된 적분법입니다.

("A method of computation for structural dynamics", Journal of the Engineering Mechanics Division, 85 (EM3): 67–94)

 

implicit method란 현재와 미래 시간의 시스템의 상태로부터

 

미래 시간의 상태를 계산하는 방법입니다.

 

반대말은 explicit method로, 현재 시간의 시스템의 상태로부터

 

미래 시간의 상태를 계산하는 방법이 있습니다.

 

 Newmark method 유도

일반적으로 운동방정식은 2차 미분방정식의 형태로 나타낼 수 있으며,

 

식으로 나타내면 아래와 같습니다.

 

$Mu''+Cu'+Ku=f(t)$

 

여기에서부터 시작해서 우리는 $u''$과 $u'$을 구하여 대입할겁니다!

 

$u(t+△t)$와 $u'(t+△t)$의 Taylor's series를 사용할 건데요.

 

Taylor's series가 무엇인지는 아래 포스팅에 설명해두었습니다.

[공업수학] Taylor's series(테일러 급수)

 

 

Newmark method는 3th order항까지 사용합니다. 그 이상의 고차항은 error term이 되죠.

 

$u_{n+1}=u_n+u'_n\triangle t+u''_n\frac{(\triangle t)^2}{2!}+u'''_n\frac{(\triangle t)^3}{3!}+...$

 

$u'_{n+1}=u'_n+u''_n\triangle t+u'''_n\frac{(\triangle t)^2}{2!}+...$

 

 

여기서 4차항 이상은 절삭하고, 3차항의 계수를 각각 $\beta$, $\gamma$로 둡니다. 즉,

 

$u_{n+1}=u_n+u'_n\triangle t+u''_n\frac{(\triangle t)^2}{2!}+\beta u'''_n (\triangle t)^3$

 

$u'_{n+1}=u'_n+u''_n\triangle t+\gamma u'''_n (\triangle t)^2$

 

후에 파라미터에 대해 설명하겠지만,

 

테일러 급수에서 나온 값이므로 $\beta$와 $\gamma$는 0에서 1사이여야 합니다.

 

여기서 함수가 선형이라면 $u'''_n$은 아래와 같이 나타낼 수 있습니다.

 

$u'''_{n+1}=(u''_{n+1}-u''_n)/\triangle t\tag{1}$

 

위 식을 Taylor's series에 의해 도출된 식에 각각 대입하여 정리하면 아래와 같습니다.

 

$u_{n+1}=u_n+u'_n\triangle t+(1-2\beta)\frac{(\triangle t)^2}{2}u''_n+\beta (\triangle t)^2 u''_{n+1}\tag{2}$

 

$u'_{n+1}=u'_n+\triangle t(1-\gamma)u''_n+\triangle t\gamma u''_{n+1}\tag{3}$

 

이 때, $0 \leq 2\beta \leq 1$, $0 \leq \gamma \leq 1$ 입니다.

 

이제 $u''$과 $u'$을 구했으니 원래의 운동방정식으로 돌아가봅시다.

 

$Mu''_{n+1}+Cu'_{n+1}+Ku_{n+1}=f(t)$

 

→ $u''_{n+1}=M^{-1}(f(t)-Cu'_{n+1}-Ku_{n+1})$

 

식 (2)을 이용해 $u''_{n+1}$을 구해봅시다.

 

$\beta (\triangle t)^2 u''_{n+1} = u_{n+1} - u_n - u'_n\triangle t - (1-2\beta)\frac{(\triangle t)^2}{2}u''_n$

 

→  $u''_{n+1} = (u_{n+1}-u_n)/{\beta (\triangle t)^2} - u'_n/{\beta\triangle t} - (1/(2\beta)-1) u''_n\tag{4}$

 

여기서 구한 $u''_{n+1}$을 식 (3)에 대입하여 $u'_{n+1}$도 구해줄 겁니다.

 

→ $u'_{n+1} = \gamma (u_{n+1}-u_n)/(\beta \triangle t) + u'_n(1-\gamma/\beta) + \triangle t u''_n(1-\gamma/(2\beta))\tag{5}$

 

그럼 운동방정식에 필요한 $u''$, $u'$, $u$가 모두 있으니 대입해볼까요?

 

$M[(u_{n+1}-u_n)/{\beta (\triangle t)^2} - u'_n/{\beta\triangle t} - (1/(2\beta)-1) u''_n]$
$+C[\gamma (u_{n+1}-u_n)/(\beta \triangle t) + u'_n(1-\gamma/\beta) + \triangle t u''_n(1-\gamma/(2\beta))]$
$+Ku_{n+1} = f(t_{n+1})\tag{6}$

 

식 (6)을 아래와 같이 explicit equation의 형태로 나타내면 어떨까요?

 

$Au_{n+1}=B_n$

 

→$u_{n+1}=A^{-1}B_n\tag{7}$

 

이 때, $A$와 $B_n$은 아래와 같습니다.

 

$A=\frac{M}{\beta(\triangle t)^2}+\frac{\gamma C}{\beta (\triangle t)}+K$

 

$B_n=f(t_{n+1})+M[\frac{u_n}{\beta(\triangle t)^2}+\frac{u'_n}{\beta(\triangle t)}+(\frac{1}{2\beta}-1)u''_n]$

       $+C[\frac{\gamma u_n}{\beta(\triangle t)}-u'_n(1-\frac{\gamma}{\beta})-\triangle t u''_n(1-\frac{\gamma}{2\beta})]$

 

 

따라서 우리는 식 (7)과 같은 방법으로 n번째 time step의 상태를 이용해

 

n+1번째 time step의 상태를 계산할 수 있습니다.