Skip to content

13. 수치해석(Numerical Analysis)

Jump to bottom
choiyeonseok edited this page Nov 20, 2019 · 1 revision

1. 이분법

주어진 범위[lo, hi] 내에서 어떤 함수 f(x)의 값이 0이 되는 지점을 수치적으로 찾아내는 기법이다.
이분법을 사용하기 위해서는 그래프상에서 x축 윗부분에서의 한 점과 아랫부분에서의 한점을 찾아야 한다.

구간의 크기를 매번 절반으로 줄여가며 정답(해)이 포함되어 있는지 확인한다.

def f(x):
    return x^2 - 4

def bisection(lo, hi):
    if f(lo) > 0:
        lo, hi = hi, lo
        
    while fabs(hi - lo) > 0.00000002:
        mid = (lo+hi) // 2
        fmid = f(mid)
        
        if fmid <= 0:
            lo = mid
        else:
            hi = mid
            
    return (lo + hi) // 2

◈ 절대 오차를 이용한 종료 판정

  • hi와 lo의 절대 오차를 이용한 반복문의 종료 판정을 구현했다. hi와 lo가 충분히 가까워질 경우 [lo, hi]구간 내에서 임의의 답을 선택해서 반환하는 것이다.

◈ 상대 오차를 이용한 종료 판정 

def infinite_bisection():
    lo = 123456123456.1234588623046875 #double
    hi = 123456123456.1234741210937500 #double
    
    while fabs(hi - lo) > 0.0000002:
        hi = (lo + hi) // 2.0
    print("finished")
  • 위 방식에는 문제가 있다. 숫자의 절대 값이 커지면 double변수에 표현할 수가 없기 때문에 절대 오차 hi값과 같게 되고 hi는 항상 2*10^-7보다 크기 때문에 무한루프가 발생하게 된다.
  • 따라서 프로그래밍 대회에서는 상대오차를 허용하고 있는데, 반환 값의 오차가 10^-7보다 크더라도, 정답 A에 대해 [Ax(1-10^-7), Ax(1+10^-7)]범위에 포함되는 답을 정답으로 인정한다.
  • 하지만 이 상대오차를 이용하는 방법보다 더 좋은 방법이 있다...!

◈ 정해진 횟수만큼 반복하기

while문을 적당한 for문으로 대체해서 반복문이 항상 정해진 만큼만 실행되도록 하는 것이다. 만약 반복문을 100번 수행하게 되면 절대 오차는 |lo-hi|/2^101이 되고 2^101은 대략 서른 한자리 수로 |hi-lo|가 대략 10^20미만의 수라면 이 오차는 항상 10^-7보다 작다.

◈ 이분법 예제1 : 단변수 다항 방정식의 수치적 해법

  • 다항방정식을 푸는 열쇠는 다항식의 두 극점 사이에는 최대 하나의 근만이 있다는 사실을 깨닫는 것이다.
  • 극점은 기울기의 부호가 바뀌는 지점이므로, 극점 사이에는 방정식의 값이 순증가하거나 순감소해야 하기 때문이다.
  • 따라서 '극점'의 위치를 알아 낼 수 있으면 이분법을 사용하여 근을 찾을 수 있다.

즉, 주어진 다항식을 미분한 뒤 재귀 호출로 해당 다항식이 0이 되는 지점을 찾고, 여기서 얻은 극점들 사이에서 이분법을 이요해 방적식의 근을 찾는 알고리즘을 작성할 수 있다.

def differentiate(poly):
    # 단일 변수 다항식 p의 계수가 주어질 때 미분한 결과 p`의 계수를 반환한다.

def solveNaive(poly):
    # 1차 혹은 2차 방정식을 푼다.

def evaluate(poly, x0):
    # 다항식의 f(x)의 계수 poly가 주어질 때, f(x0)를 반환한다.
    
# 방정식 해의 절대 값은 L이하여야 한다.
L = 25

# 방정식 sum(poly[i] * x^i) = 0의 근을 반환한다.
def solve(poly):
    
    n = len(poly) - 1
    
    # 기저 사례로서 2차 이하의 다항식은 그냥 풀수 있다.
    if n <= 2:
        return solveNaive(poly)
    
    # 방정식을 미분해서 n-1차 방정식을 얻는다.
    derivative = differentiate(poly)
    sols = solve(derivative)
    
    # 극점들 사이를 하나하나를 검사하며 근이 있나 확인한다.
    sols.insert(0, -L-1)
    sols.insert(len(sols), L+1)
    
    ret = 0
    for i in range(len(sols) - 1):
        x1, x2 = sols[i], sols[i+1]
        y1, y2 = evaluate(poly, x1), evaluate(poly, x2)
        
        # f(x1)과 f(x2)의 부호가 같은 경우 답은 없다.
        if y1*y2 > 0:
            continue
        
        # 불변조건 : f(x1) <= 0 < f(x2)
        if y1 > y2:
            y1, y2 = y2, y1
        
        #이분법 사용
        for iter in range(100):
            mx = (x1 + x2) / 2
            my = evaluate(poly, mx)
            if y1*my > 0:
                y1, x1 = my, mx
            else:
                y2, x2 = my, mx
        ret.append((x1 + x2) / 2)
    ret.sort()
    return ret

◈ 이분법 예제2 : 전세금 균등상환

  • 전세금 N원을 연 연이율 P%로 대출받고, 이를 M(1<=M<=120)개월 동안 매달 일정액 C원씩 갚는다고 설정하자.

  • 대출의 잔금은 대출 금액 N원에서 시작한다.

  • 한 달이 지날 때마다 대출 잔금이 월 이자 P/12%만큼 불어난다.

  • 이자가 추가 된 다음 월 상환액 C를 대출 잔금에서 제한다.

  • 이때, M개월에 걸려 모든 대출 금액을 다 갚기 위해서는 한달에 최소 얼마씩 갚아야 할까?

매달 C원씩 대출을 상환할 때, M개월 후의 대출 잔액을 Balance(C)라고 한다면 이 문제는 단순히 balance(c)= 0의 방정식을 푸는 문제가 된다.

# 전세금 균등상환 문제를 해결하는 이분법의 구현

# amount원을 연 이율 rates퍼센트로 duration 월 간 한달에 monthly Payment로 남았을 때 대출 잔금은??
def balance(amount, duration, rates, monthlyPayment):
    balance = amount
    for i in range(duration):
        # 이자가 붙는다.
        balance *= 1.0 + (rates/12.0) / 100.0
        
        # 상환액을 잔금에서 제한다.
        balance -= monthlyPayment
        
    return balance

# amount 원을 연 이율 rates퍼센트로 duration월 간 갚으려면 한달에 얼마씩 갚아야하나?
def payment(amount, duration, rates):
    
    lo = 0  # 돈을 아예 갚지 않는 경우    
    hi = amount*(1.0 + rates/12.0) / 100.0 # 한달 후에 전액을 상환하는 경우 
    for iter in range(100):
        if balance(amount, duration, rates, mid) <= 0:
            hi = mid
        else:
            lo = mid
    return hi

Clone this wiki locally