[Python] 13705번 Ax+Bsin(x)=C

https://www.acmicpc.net/problem/13705

13705번: Ax+Bsin(x)=C

---

 

24/02/14

 

 

아주 유명한 문제다. 사용하는 알고리즘은 매우 단순하다.

 

그냥 이분 탐색에 약간의 수학적 지식.

 

이 문제가 그렇게 악명이 높은 이유는 매우 높은 정밀도를 요구하기 때문이다.

 

파이썬을 사용하면 쉽게 구현할 수 있다.

---

 

문제 접근 방식:

 

 

함수 $f(x) = Ax + B\mathrm{sin}(x) - C$이 $B \leq A$조건 하에서 단조 증가함수임을 보이기만 한다면, 일단 이분 탐색을 사용할 수 있는 조건이 성립된다.

 

함수 $f$를 미분한 결과는 $f'(x) = A + B\mathrm{cos}(x)$이다. 이때, $-1 \leq \mathrm{cos}(x) \leq 1$이 성립하므로, $B \leq A$조건 하에서는 $f'(x) \leq 0$이다.

 

따라서, 단조 증가 함수이므로 근을 가진다면 오직 하나만 가짐을 보장할 수 있다.

 

근을 찾는 방법은 이분 탐색을 통해 구현할 수 있다.

 

만약 $f(\textrm{mid}) > 0$이라면, 근은 $\textrm{mid}$값보다 왼쪽에 존재함을 알 수 있다.

 

따라서 이 경우에는 $\textrm{end}$를 $\textrm{mid}$로 다시 설정하여 이분 탐색을 진행하면 된다.

 

만약 $f(\textrm{mid}) < 0$이라면, 근은 $\textrm{mid}$값보다 오른쪽에 존재하므로 $\textrm{start}$를 $\textrm{mid}$로 설정하여 이분 탐색을 진행하면 된다.

 

두 번째로 필요한 지식은 사인 함수의 테일러 전개(정확히 말하면 맥클로린 전개) 형태이다.

 

사인 함수의 테일러 전개는 다음과 같은 형태를 가지고 있다.

 

$$\mathrm{sin}(x) = x - \frac{x^3}{3!} + \frac{x^5}{5!} - \cdots$$

 

이를 그대로 코드로 구현하면 된다.

 

파이썬의 decimal모듈이 필수적인데, 공식 문서를 보면 sin함수에 대한 구현 예시도 나와 있어서, 직접 구현하기 귀찮다면 복붙해도 괜찮다.

 

https://docs.python.org/3/library/decimal.html

decimal — Decimal fixed point and floating point arithmetic

 

정밀도는 처음에는 100자리로 했으나, 시간초과가 자꾸 떠서 50자리로 줄였다.

 

이후 중간에 무한 루프가 날 수 있는 케이스가 있을 수 있겠다 싶어서 무한 루프 방지용 코드를 삽입해주고, 반올림도 잘해주니 맞았습니다를 받을 수 있었다.

 

decimal모듈에 대한 설명은 위의 공식 문서를 참조하도록 하자. 

---

아래는 내가 위의 접근 방식과 같이 작성한 파이썬 코드이다. 더보기를 누르면 확인할 수 있다.

# 13705번 Ax+Bsin(x)=C
# 수치해석
import sys
input = sys.stdin.readline
from decimal import *
getcontext().prec = 50
getcontext().rounding=ROUND_HALF_UP
pi = Decimal('3.1415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421170679')
err = Decimal('0.0000000000000000000000000000000000000000000001')
def sin(x):
    x = x % (Decimal('2')*pi)
    i, lasts, s, fact, num, sign = 1, 0, x, 1, x, 1
    no_inf_loop = 0
    while abs(s-lasts) > err and no_inf_loop < 1000:
        lasts = s
        i += 2
        fact *= i * (i-1)
        num *= x * x
        sign *= -1
        s += num / fact * sign
        no_inf_loop += 1
    return +s
a, b, c = map(Decimal, input().rstrip().split())
start, end = Decimal('0'), Decimal('150000')
no_inf_loop = 0
while start < end and (end-start) > err and no_inf_loop < 1000:
    mid = (start + end) / Decimal('2')
    ans = a*mid + b*sin(mid) - c
    if ans > 0:
        end = mid
    elif ans < 0:
        start = mid
    else:
        break
    no_inf_loop += 1

print(f'{round(mid, 6):.6f}')