[백준] 11401: 이항 계수 3 - 분할 정복, 페르마의 소정리

다른 코드 참조: O, 푼 횟수: 8

들어가기에 앞서..

이항 계수를 1,000,000,007로 나눈 나머지를 구하는 문제로  ₙCₖ % 1,000,000,007을 출력하면 되는 간단한 문제 같지만 다음과 같은 경우에 에러가 발생한다.

 

1. 반복문으로 푸는 경우

먼저 반복문으로 푸는 경우를 생각해 볼 수 있다.

import sys

n, k = map(int, sys.stdin.readline().split())
fac = [1 for i in range(n+1)]

for i in range(1, n+1):
	fac[i] = fac[i-1] * i

print((fac[n]//(fac[k]*fac[n-k])) % 1000000007)

문제는 n≤4,000,000 이라는 점이다.

n이 20만 돼도 20! = 51,090,942,171,709,440,000 이므로 숫자가 너무 커서 속도가 느려진다.

 

그러면 팩토리얼을 구하는 과정에서 계속해서 나머지 연산을 하여 1,000,000,007의 나머지가 저장되도록해서 크기가 너무 커지게 하지 않는 경우가 있다.

import sys

n, k = map(int, sys.stdin.readline().split())
fac = [1 for i in range(n+1)]

for i in range(1, n+1):
	fac[i] = (fac[i-1] * i) % 1000000007

print((fac[n]//(fac[k]*fac[n-k])) % 1000000007)

 

하지만 나눗셈은 모듈러 연산의 분배법칙이 적용되지 않는다. 따라서 p = 1,000,000,007 일 때

 

2. 재귀로 푸는 경우

재귀로 푸는 경우에도 n이 너무 크기 때문에 스택 오버플로우가 발생한다.

 

따라서 시간초과와 메모리 초과를 일으키지 않고 문제를 풀려면, 팩토리얼을 구하는 과정에서 계속해서 나머지 연산을 하여 1,000,000,007의 나머지가 저장되도록해서 크기가 너무 커지게 하지 않아야 한다. 하지만 1번의 경우에서 봤듯이 모듈러 연산에서 나눗셈이 없다는 점이다. 따라서 나눗셈을 곱셈으로 바꿔서 계산할 수 있는 페르마의 소정리가 필요하다.

 

페르마의 소정리

p가 소수이면, 모든 정수 a에 대해  aᴾ≡ a (mod p)이다.

혹은

p가 소수이고 a가 p의 배수가 아니면,  aᴾ⁻¹ ≡ 1 (mod p)이다.

 

증명

1. 기약잉여계를 사용한 증명

1. a와 p가 서로소라고 하자.
2. A = {1*a, 2*a, 3*a,... , (p-1)*a}인 집합이 있다고 가정할 때, 당연하게도 각 원소들의 mod p는 서로 다르다.
3. 또한, 1번에 의해 p의 배수인 원소가 없으므로 A의 원소들을 p로 나눈 나머지들의 집합을 B라고 했을 때, B = {1, 2,..., p-1}
4. 따라서 (A의 각 원소들의 mod p의 곱) ≡ (B의 각 원소들의 곱) (mod p)
5. (1*a mod p)*(2*a mod p)*...*((p-1)*a mod p) ≡ 1*2*...*(p-1) (mod p)
6. (1*a)*(2*a)*...*((p-1)*a) ≡ 1*2*...*(p-1) (mod p)
7. aᵖ⁻¹*(p-1)! ≡ (p-1)! (mod p)
8. p는 소수이므로 (p-1)! 과 p는 서로소이다. 따라서 양변에서 (p-1)!를 나누어 줄 수 있다. ((p-1)! 과 p가 서로소가 아니라면 mod p가 0이 될 수 있으므로 나눠 줄 수 없다)
∴ aᴾ⁻¹ ≡ 1 (mod p)

aᵖ ≡ a (mod p)는 수학적 귀납법을 통해 증명할 수 있다.

 

aᴾ⁻¹ ≡ 1 (mod p)라는 것을 알았으므로 다음과 같이 유도할 수 있다. aᴾ⁻² ≡ a⁻¹ (mod p)

따라서, 분수꼴을 곱셈꼴로 바꿀 수 있다.

 

풀이과정

앞에서 정리한 것들을 종합하면 문제를 풀 수 있다. 

n! * (k! * (n-k)!)⁻¹ ≡ n! * (k! * (n-k)!)ᴾ⁻² (mod p)

곱셈으로 변환했으므로 각각 나눠서 나머지를 구한 후 곱해도 무방하다!

(n! % p) * ((k! % p) * ((n - k)! % p))ᴾ⁻²  % p = n! * (k! * (n - k)!)⁻¹ % p

 

코드

import sys


n, k = map(int, sys.stdin.readline().split())
fac = [1 for i in range(n+1)]
p = 1000000007


def power(a, x):
    if x == 1:
        return a % p

    r = power(a, x//2)
    if x % 2 == 0:
        return (r * r) % p
    else:
        return (r * r * a) % p


for i in range(1, n+1):
    fac[i] = (fac[i-1]*i) % p		#계속해서 나머지를 저장해 숫자가 커지지 않게 한다.

bot = (fac[k]*fac[n-k]) % p
bot = power(bot, p-2)			#(k!*(n-k)!)^(p-2)를 bot에 저장한다.
print((fac[n]*bot) % p)