logo

뮬러 메소드 📂수치해석

뮬러 메소드

메소드

f(α)=0f (\alpha) = 0 이라고 하자. 초기값 x0,x1,x2x_{0} , x_{1} , x_{2} 과 wn:=f[xn,xn−1]+f[xn,xn−2]−f[xn−2,xn−1] w_{n} := f [x_{n} , x_{n-1} ] + f [ x_{n} , x_{n-2} ] - f [ x_{n-2} , x_{n-1} ] 에 대해 xn+1:=xn−2f(xn)wn±wn2−4f(xn)f[xn,xn−1,xn−2] x_{n+1} : = x_{n} - {{ 2 f ( x_{n} ) } \over { w_{n} \pm \sqrt{ w_{n}^{2} - 4 f (x_{n} ) f [ x_{n} , x_{n-1} , x_{n-2} ] } }} 과 같이 정의된 수열 {xn}\left\{ x_{n} \right\} 은 n→∞n \to \infty 일 때 α\alpha 로 p≈1.84p \approx 1.84 차 수렴한다.

  • 단, (wn±wn2−4f(xn)f[xn,xn−1,xn−2])∈C\left( w_{n} \pm \sqrt{ w_{n}^{2} - 4 f (x_{n} ) f [ x_{n} , x_{n-1} , x_{n-2} ] } \right) \in \mathbb{C} 는 ++ 와 −- 둘 중 ∣wn±wn2−4f(xn)f[xn,xn−1,xn−2]∣\left| w_{n} \pm \sqrt{ w_{n}^{2} - 4 f (x_{n} ) f [ x_{n} , x_{n-1} , x_{n-2} ] } \right| 이 큰 쪽으로 선택한다.
  • 복소수의 절댓값은 ∣x+iy∣:=x2+y2|x + iy| := \sqrt{ x^2 + y^2} 와 같이 계산된다.

설명 1

뮬러 메소드는 시컨트 메소드의 응용으로 볼 수 있는데, 여전히 뉴턴-랩슨 메소드보다 느리지만 복소수해를 찾을 수 있다는 특징이 있다. 미분이나 부호를 이용한 풀이가 아니라 근의 공식을 이용해서 풀기 때문에 기하학적인 의미에 매달릴 필요가 없으며 복소수를 커버한다.

예시

주어진 (xi,f(xi))i=n−1n+1\left( x_{i} , f (x_{i}) \right)_{i = n-1}^{n+1} 에 대해 쿼드러틱 폴리노미얼 p2(x)=f(xn+1)+(x−xn+1)f[xn+1,xn]+(x−xn+1)(x−xn)f[xn+1,xn,xn−1] p_{2} (x) = f(x_{n+1} ) + (x - x_{n+1} ) f [ x_{n+1} , x_{n} ] + ( x - x_{n+1} ) ( x - x_{n} ) f [ x_{n+1} , x_{n} , x_{n-1} ] 은 유일하게 존재하며 p(xi)=f(xi)p ( x_{i} ) = f(x_{i} ) 을 잘 만족한다. 이를 (x−xn+1)(x - x_{n+1}) 에 대한 이차식으로 바꿔보자.

f[xn+1,xn,xn−1]=f[xn+1,xn−1]−f[xn−1,xn]xn+1−xn  ⟹  (xn+1−xn)f[xn+1,xn,xn−1]=f[xn+1,xn−1]−f[xn−1,xn]  ⟹  (x−xn)f[xn+1,xn,xn−1]=f[xn+1,xn−1]−f[xn−1,xn]+(x−xn+1)f[xn+1,xn,xn−1]  ⟹  f[xn+1,xn]+(x−xn)f[xn+1,xn,xn−1]=f[xn+1,xn]+f[xn+1,xn−1]−f[xn−1,xn]+(x−xn+1)f[xn+1,xn,xn−1]  ⟹  f[xn+1,xn]+(x−xn)f[xn+1,xn,xn−1]=wn+(x−xn+1)f[xn+1,xn,xn−1]  ⟹  (x−xn+1)f[xn+1,xn]+(x−xn)(x−xn+1)f[xn+1,xn,xn−1]=wn(x−xn+1)+(x−xn+1)2f[xn+1,xn,xn−1] \begin{align*} & f [ x_{n+1} , x_{n} , x_{n-1} ] = { { f [x_{n+1} , x_{n-1} ] - f [ x_{n-1} , x_{n} ] } \over { x_{n+1} - x_{n} }} \\ \implies& ( x_{n+1} - x_{n} ) f [ x_{n+1} , x_{n} , x_{n-1} ] \\ & = f [x_{n+1} , x_{n-1} ] - f [ x_{n-1} , x_{n} ] \\ \implies& ( x - x_{n} ) f [ x_{n+1} , x_{n} , x_{n-1} ] \\ & = f [x_{n+1} , x_{n-1} ] - f [ x_{n-1} , x_{n} ] + ( x - x_{n+1} ) f [ x_{n+1} , x_{n} , x_{n-1} ] \\ \implies& f [x_{n+1} , x_{n} ] + ( x - x_{n} ) f [ x_{n+1} , x_{n} , x_{n-1} ] \\ & = f [x_{n+1} , x_{n} ] + f [x_{n+1} , x_{n-1} ] - f [ x_{n-1} , x_{n} ] + ( x - x_{n+1} ) f [ x_{n+1} , x_{n} , x_{n-1} ] \\ \implies& f [x_{n+1} , x_{n} ] + ( x - x_{n} ) f [ x_{n+1} , x_{n} , x_{n-1} ] \\ & = w_{n} + ( x - x_{n+1} ) f [ x_{n+1} , x_{n} , x_{n-1} ] \\ \implies& ( x - x_{n+1} ) f [x_{n+1} , x_{n} ] + ( x - x_{n} ) ( x - x_{n+1} ) f [ x_{n+1} , x_{n} , x_{n-1} ] \\ & = w_{n} ( x - x_{n+1} ) + ( x - x_{n+1} )^2 f [ x_{n+1} , x_{n} , x_{n-1} ] \end{align*}

따라서 p2(x)=f(xn+1)+(x−xn+1)f[xn+1,xn]+(x−xn+1)(x−xn)f[xn+1,xn,xn−1]=f(xn+1)+wn(x−xn+1)+(x−xn+1)2f[xn+1,xn,xn−1] \begin{align*} p_{2} (x) =& f(x_{n+1} ) + (x - x_{n+1} ) f [ x_{n+1} , x_{n} ] + ( x - x_{n+1} ) ( x - x_{n} ) f [ x_{n+1} , x_{n} , x_{n-1} ] \\ =& f(x_{n+1} ) + w_{n} (x - x_{n+1} ) + ( x - x_{n+1} )^2 f [ x_{n+1} , x_{n} , x_{n-1} ] \end{align*}

x=xn+1x = x_{n+1} 을 대입하면 p2(xn+1)=f(xn+1)+0+0p_{2} (x_{n+1}) = f(x_{n+1}) + 0 + 0 이 되므로, p2(x)=0p_{2} (x) = 0 이 되도록 하는 (x−xn+1)(x - x_{n+1}) 을 찾으면 결국 f(xn+1)=0f(x_{n+1} ) = 0 을 찾는 것이나 마찬가지다. 원래 구하고 싶은 것이 lim⁡n→∞f(xn)=0\displaystyle \lim_{n \to \infty} f ( x_{n} ) = 0 를 만족하는 xnx_{n} 이었으므로 p2(x)=f(xn+1)+wn(x−xn+1)+(x−xn+1)2f[xn+1,xn,xn−1] p_{2} (x) = f(x_{n+1} ) + w_{n} (x - x_{n+1} ) + ( x - x_{n+1} )^2 f [ x_{n+1} , x_{n} , x_{n-1} ] 의 (x−x2)(x - x_{2}) 에 대해 근의 공식을 사용하면 된다. 근의 공식으로 xn+1−xn=−wn±wn2−4f(xn)f[xn,xn−1,xn−2]2f[x2,x1,x0] x_{n+1} - x_{n} = {{ - w_{n} \pm \sqrt{ w_{n}^{2} - 4 f (x_{n} ) f [ x_{n} , x_{n-1} , x_{n-2} ] } } \over { 2 f [ x_{2} , x_{1} , x_{0} ] }} 를 얻는데, 실제 계산을 컴퓨터가 하기 때문에 유효숫자 등의 문제 때문에 유리화를 거꾸로 해줘서 xn+1=xn−2f(xn)wn±wn2−4f(xn)f[xn,xn−1,xn−2] x_{n+1} = x_{n} - {{ 2 f ( x_{n} ) } \over { w_{n} \pm \sqrt{ w_{n}^{2} - 4 f (x_{n} ) f [ x_{n} , x_{n-1} , x_{n-2} ] } }} 을 공식으로 사용한다.

한편 f′f ' 가 꼭 존재할 필요는 없으나, f’(α)≠0f’ ( \alpha ) \ne 0 이라면 n→∞n \to \infty 일 때 wn→f′(α)w_{n} \to f ' ( \alpha ) 이기 때문에 미분계수를 찾는 용도로도 쓸 수 있기는 하다.

구현

다음은 뮬러 메소드를 R 로 구현한 것이다. 실제 정리에선 세 점 (x0,x1,x2)(x_{0} , x_{1} , x_{2}) 을 가지고 시작하지만 실제론 x0,x1x_{0} , x_{1} 만 있어도 x2:=x0+x12\displaystyle x_{2} := {{x_{0} + x_{1}} \over {2}} 을 결정할 수 있어서 두 개의 초기값만 받아서 사용한다. 이 때 (x0,x1)(x_{0} , x_{1}) 이 실제로 실수일지라도 자료형은 복소수여야한다.

dd<-function(f,X){
  if(length(X)==1) {return(f(X))}
  temp<-numeric(0)
  for(i in 1:length(X)) {temp[i]<-prod(X[i]-X[-i])}
  return(sum(f(X)/temp))
}
 
Muller<-function(f,x0,x1)
{
  x2<-(x0+x1)/2
  
  while(Mod(f(x2))>10^(-3)){
    w = dd(f,c(x2,x1)) + dd(f,c(x2,x0)) - dd(f,c(x0,x1))
    p = w + sqrt(w^2 - 4*f(x2)*dd(f,c(x2,x1,x0)))
    m = w - sqrt(w^2 - 4*f(x2)*dd(f,c(x2,x1,x0)))
    
    if(Mod(p)>Mod(m)) {w<-p} else {w<-m}
    
    x0<-x1
    x1<-x2
    x2 <- (x2 - (( 2*f(x2) ) / w))
  }
  return(x2)
}
 
f<-function(x) {x^2 + 9}
Muller(f,-6+0i,-5+0i)
 
g<-function(x) {x^2 + x + 1}
Muller(g,-1i,-2+0i)

다음은 위 코드를 실행 시킨 결과다. f(x)=x2+9=f(x) = x^2 + 9 = 의 허근 중 하나인 x=−3ix = -3 i 와 g(x)=x2+x+1=0g(x) = x^2 + x + 1 = 0 의 복소근 중 하나인 x=−1+3i2\displaystyle x = - {{1 + \sqrt{3} i } \over {2}} 을 잘 찾아냈음을 확인할 수 있다.

20190321\_165937.png

  1. Atkinson. (1989). An Introduction to Numerical Analysis(2nd Edition): p73~74. ↩︎