문제


스포 방지선

___________________________________________________________________________________________________________












___________________________________________________________________________________________________________


풀이


Euclidean Algorithm의 느낌으로 생각하면, \(ax+by=s\)를 만족하는 음이 아니고 서로소인 두 정수 \(x, y\)를 찾으면 된다.

(UPD: 사실 이 점을 증명하는 것도 trivial 한 것은 아니다. 머리를 좀 쓰긴 해야함) 

일반적으로 \(ax+by=s\)를 푸는 것은 Extended Euclidean Algorithm을 알면 쉽게 구현 가능하다. 

\(\text{gcd}(a, b)=1\)이 되도록 식을 reduce 시키고, 얻은 새로운 식을 \(ax+by=s\)라고 재정의하자. 

이제 \(ax+by=1\)을 만족하는 정수 \(x, y\)를 Extended Euclidean Algorithm (Modular Inverse)으로 구할 수 있다. 

위 식을 단순하게 \(s\)배 하기만 하면, \(ax+by=s\)인 정수 \(x, y\)를 하나 찾을 수 있다. 이를 \(x_0, y_0\)라 하자.


우리는 \(ax + by=s\)가 성립하면, 정수 \(k\)에 대해 \(x = x_0 + bk\), \(y = y_0 - ak\)라고 쓸 수 있음을 안다. 

이제 \(x_0, y_0\)를 적당히 이 방식대로 바꿔서, \(x_0 < b\)가 성립하도록 할 수 있다. 

다시 \(x = x_0 + bk\), \(y= y_0-ak\)로 solution set을 parametrize하자.

\(x, y \ge 0\)이어야 하니, \(k \ge 0\)인 경우만 보면 된다는 것을 알 수 있다. 


이제 목표는 \(k \ge 0\), \(x_0+bk \ge 0\), \(y_0-ak \ge 0\), \(\text{gcd}(x_0+bk, y_0-ak)=1\)인 정수 \(k\)가 존재하는지 여부다. 

꽤 어려운 문제처럼 보이지만, 알고리즘 문제를 푸는 입장에서 solution은 간단하다. 

단순히 \(k\)에 대해서 iterate 하면서, 조건을 만족하는 \(k\)가 나오면 YES를 찍고, 아니면 NO를 찍으면 된다. 

\(k \ge 0\), \(x_0+bk \ge 0\), \(y_0-ak \ge 0\)를 만족하는 정수 \(k\)의 개수는 \(10^{18}\)-scale까지 커질 수 있다는 것을 생각하자.

이 알고리즘이 시간 안에 답을 낸다는 사실은 (구현하면 0ms 뜬다) 수학적으로 전혀 자명하지 않다.


1년 반 전에 이 문제를 풀었을 때 이 점은 개인적으로 숙제로 남아 있었다.

오늘 이 풀이가 시간 안에 도는 이유를 아는 분에게 질문 받아서, 고민을 하는 시간을 가져 해답을 찾았다.

사실 이 블로그를 쓰는 이유도 풀이를 쓰는 것보단 이 해답을 기록하기 위한 것이다. 


Claim: 문제 조건에서 \(\text{gcd}(x_0+bk, y_0-ak)=1\)을 만족하는 최소의 음이 아닌 정수 \(k\)는 최대 \(2^{15}\)이다. 


이 Claim이 증명되면, 풀이가 시간 안에 도는 이유는 자명하다. 

\(k\)에 대한 iteration이 최대 \(2^{15}\)번 돌기 때문에, 코드의 시간복잡도는 대강 \(\mathcal{O}(2^{15} \log s)\) 정도가 된다.

이 정도의 시간복잡도/상수면, 코드가 0ms 안에 도는 것도 충분히 설명된다. 


Claim 증명


자연수 \(D\)가 있어, 각 \(k=0, 1, \cdots D-1\)에 대해 \(\text{gcd}(x_0+bk, y_0-ak) \neq 1\)이라 하자.

최대공약수가 \(1\)이 아니라는 것은 공통 소인수가 있다는 것이다. 

그러니 각 \(k=0, 1, \cdots , D-1\)에 대해 \(x_0+bk\)와 \(y_0-ak\)의 공통 소인수를 아무거나 하나 잡아서 이를 \(p_k\)라 하자.  


우선 \(s \equiv ax_0 + by_0 \equiv a(x_0+bk)+b(y_0-ak) \equiv 0 \pmod{p_k}\)이므로, \(p_k|s\)를 얻는다.

즉, \(p_k\)로 가능한 소수의 개수는 유한하다. (특히, \(s \le 10^{18}\)이므로, 최대 15개가 존재한다)


\(P = p_i = p_j\)인 두 정수 \(0 \le i, j < D\)가 있다고 가정하자.

그러면 \(x_0 + bi \equiv x_0 + bj \equiv y_0 - ai \equiv y_0 - aj \equiv 0 \pmod{P}\)가 성립한다.

이를 정리하면, \(P|b(i-j)\), \(P|a(i-j)\)가 성립한다. \(\text{gcd}(a,b)=1\)이므로, \(i \equiv j \pmod{P}\)를 얻는다. 


이 시점에서 상황을 정리해보자. 

  • \(s\)의 소인수는 최대 15개가 있다. 
  • 각 \(k=0, 1, \cdots D-1\)에 대해, 소수 \(p_k\)는 \(s\)의 소인수다.
  • \(s\)의 각 소인수 \(P\)에 대해서, \(P = p_k\)인 index \(k\)의 수열은 공차가 \(P\)인 등차수열에 완전히 포함된다. 

이 상황을 약간 다르게 해석해보자.

  • \(s\)의 각 소인수 \(P\)에 대해서, \(P = p_k\)인 index \(k\)의 수열은 공차가 \(P\)인 등차수열에 완전히 포함된다. 
  • \(s\)의 각 소인수 \(P\)에 대해서, \(P = p_k\)인 index \(k\)의 수열을 포함하는 등차수열 \(\{v_P + Pn\}_{n \in \mathbb{Z}}\)를 생각하자.
  • 각 소인수에 대해서 얻은 등차수열의 합집합을 생각하면, \(0, 1, \cdots , D-1\)을 전부 cover 한다.
  • 즉, \(\{0, 1, 2, \cdots , D-1\} \subseteq \cup_{P|s, P \text{      is prime}} \{v_P + Pn\}_{n \in \mathbb{Z}} \)이다.

그런데, 중국인의 나머지 정리의 느낌으로 생각하면, \(\cup_{P|s, P \text{      is prime}} \{v_P + Pn\}_{n \in \mathbb{Z}} \)는 절대로 정수 집합 전체가 될 수 없다. 


이제 문제를 다음과 같은 느낌으로 바꾸어서 생각할 수 있다. 

  • 정수 전체를 cover 하지는 않는 등차수열이 \(k\)개가 있다. 
  • 이 등차수열들은 최대 몇 개의 연속한 자연수를 cover 할 수 있을까?

이 문제는 어느 정도 유명한 문제다. 개인적으로 아주 좋아하는 정리 중 하나를 소개한다.


Theorem: (Erdős, Crittenden, Vanden Eynden)

\(F\)는 \(k\)개의 등차수열 \(a_i + d_i \mathbb{Z}\)의 집합으로, (\(1 \le i \le k\)) \(1 < d_1 < d_2 \cdots < d_k\)다. 

만약 \(F\)가 연속한 \(2^k\)개의 정수를 cover 한다면, \(F\)는 정수 전체를 cover 하는 집합이다. 


Theorem 증명 (Taken From "Problems From The Book" pp. 152-153)

정수 \(x, x+1, \cdots x+2^k-1\)이 모두 \(F\)에 의해 cover 된다고 가정하자.

각 \(0 \le t <2^k\)에 대해서, 적당한 \(1 \le j \le k\)가 있어 \(x+t \equiv a_j \pmod{d_j}\)다. 

이를 다르게 해석하면, 각 \(0 \le t < 2^k\)에 대해서 \(\displaystyle \prod_{j=1}^k \left(1 - e^{\frac{2\pi i}{d_j} ( x+t-a_j)} \right) = 0\)라는 것이다. (동치조건이다)


이 식을 전개하면, \(\sum_{I \subset S} \alpha_I \cdot e^{(x+t) \beta_I} = 0\) 형태로 쓸 수 있음을 알 수 있다.

단, \(S = \{1, 2, \cdots k\}\)이며, \(\alpha_I\)와 \(\beta_I\)는 \(I, a_i, d_i\)에만 영향을 받는 값이다. 정확하게 말하자면, $$ \prod_{j=1}^k \left(1- e^{\frac{2\pi i}{d_j}(x+t-a_j)} \right) = \sum_{I \subset S} (-1)^{|I|} \cdot e^{-2\pi i \sum_{j \in I} a_j/d_j} \cdot e^{2\pi i (x+t) \sum_{j \in I} 1/d_j} $$이 성립한다. 이제 \(z_I = e^{\beta_I}\)라고 두면, 각 \(0 \le t < 2^k\)에 대해서 \(\sum_{I \subset S} \alpha_I z^{x+t}_I = 0\)이다. 


이제 \(u_n = \sum_{I \subset S} \alpha_I z^n_I\)라고 하자. \(u_n = 0\)인 것은, \(n\)이 \(F\)에 의해 cover 되는 것과 동치임을 정의상 알 수 있다.

식의 형태를 보면서 linear recurrence의 characteristic polynomial이 생각이 날 것이다. 

다항식 \(\displaystyle \prod_{I \subset S} (X-z_I) = X^{2^k} + A_{2^k-1} X^{2^k-1} + \cdots + A_1 X^1 + A_0\)를 생각하자.

그러면 각 \(I \subset S\)에 대해 \(z^{n+2^k}_I + A_{2^k-1} z^{n+2^k-1}_I + \cdots + A_1 z^{n+1}_I + A_0 z^n_I = 0\)이다.

이 식들을 적당히 선형결합하면, \(u_{n+2^k} + A_{2^k-1} u_{n+2^k-1} + \cdots + A_1 u_{n+1} + A_0 u_n = 0\)이다. 

한편, 이미 \(u_x = u_{x+1} = \cdots = u_{x+2^k-1} = 0\)을 알고 있으니, 모든 \(n\)에 대해서 \(u_n = 0\)임이 자명하다. 

증명 과정에서는 \(A_0 \neq 0\)이라는 (자명한) 사실이 사용된다. 나머지는 전형적인 귀납법. 증명 끝. 


이제 끝났다. 지금까지의 관찰을 합쳐서, \(D \le 2^{15}\)임을 알 수 있다. 증명 끝.


Note. 관련된 문제인 RMM 2010 Problem 1을 참고하면 도움이 될 것 같다. 출제자는 Dan Schwarz.



'PS > 수학 계열 PS' 카테고리의 다른 글

Project Euler 350문제  (0) 2019.12.02
Project Euler 691  (1) 2019.12.01
Project Euler 300문제  (0) 2019.11.23
9월초 수학 PS 일지  (3) 2019.09.11
7-8월의 수학 PS 일지  (3) 2019.08.28



최근에 쉬운 문제 몇 개를 더 풀어서 300문제를 찍었다. 어려운 문제 몇 개는 종강하면 포스팅을 할 듯.

'PS > 수학 계열 PS' 카테고리의 다른 글

Project Euler 691  (1) 2019.12.01
BOJ #9267 - A + B (SEERC 2013)  (2) 2019.11.24
9월초 수학 PS 일지  (3) 2019.09.11
7-8월의 수학 PS 일지  (3) 2019.08.28
Project Euler #645  (0) 2019.03.01


스코어보드


우여곡절 끝에 5등으로 첫 ICPC를 마무리했다. A Bus With No Drivers 팀으로 bjwj5505, junie, rkm0959가 한 팀이다.

예상보다는 훨씬 좋은 성적으로 대회를 끝냈지만, 여러 아쉬움도 남는 대회였다.


서울대학교가 top 6팀 중에 5개라는 것도 참 신기하다.


대회에서 우리 팀이 문제를 푼 과정을 대략 복기를 해보면,

  • 시작하자마자 bjwj5505가 J를 잡고, K를 나한테 넘겼다. 
  • J에 대한 추측을 기반으로 답안 제출 + WA
  • junie가 A가 쉽다고 해서 A를 짜서, 맞았다. 8분. 
  • junie B 읽기 시작.
  • 나는 EFG를 읽었는데 하나 같이 어려워보였다. H를 읽었는데 금광 문제와 똑같았다. 
  • 반복되는 J 뇌절
  • 그런데 알고보니 뒤에 있는 I, L번이 많이 풀려있었다. 
  • bjwj5505 L로 이동, junie B 놓고 I로 이동
  • I, L 풀이는 금방 나왔다
  • H, I, L 코딩을 동시다발적으로 해서 각각 41, 57, 51분에 맞았다. I는 1틀.
  • junie는 다시 B로, 나는 아까 읽은 K로, bjwj5505는 뇌절했던 J로 이동.
  • bjwj5505가 J에 대한 추측을 수정하면서 73분에 해결했다. bjwj5505는 D로 이동.
  • 나는 K에 대한 내 추측에 확신이 있었고, bjwj5505도 동의했다. 바로 짜서 90분에 해결했다.
  • bjwj5505가 D에 대한 추측과 알고리즘을 완성했다. 하지만 WA.
  • junie는 B 식을 짰지만 답이 잘 나오지 않았고, 내가 B로 가서 같이 DP 식과 코드를 고쳤다.
  • B를 int/ll 문제로 한 번 WA를 받고, 127분에 맞았다. 
  • D 코드를 같이 봤지만 문제가 없어보였다. 3번 WA를 받고 잠시 D에서 손을 놓기로 했다.
  • G를 같이 읽었는데, bjwj5505가 그냥 풀어버렸다. 코드 디버깅은 같이 했다. 
  • 1 WA를 받은 뒤 193분에 AC를 받았다.
  • D에 대한 확신을 많이 잃은 뒤라, F를 읽자는 의견이 나왔다.
  • 기하가 너무 빡세보여서 나는 D를 어떻게든 잡고 싶었다. stress를 돌리자는 의견이 나왔다.
  • 내가 전탐색으로 최소 막대기 개수를 구하는 코드를 짰다. (Dnaive.cpp)
  • 내가 우리 코드가 정당한 답을 내는지 확인하는 코드를 D.cpp에 추가했다.
  • junie가 테스트 제너레이터를 짰다. (Dgen.cpp)
  • junie가 stress test 코드를 짰다. (Dstress.cpp)
  • 돌렸더니 금방 반례가 나왔다. 로직의 문제도 곧 나왔다. 
  • 로직의 문제를 해결하는 코드를 냈는데, TLE를 예상했지만 바로 AC. 263분에 D를 맞았다.
  • 남은 문제는 너무 어려워보여서 던졌다. C, F가 할 만한 문제인줄은 몰랐지 ㅋㅋ
  • 끝나고 Ternion 팀하고 같이 술과 고기를 먹었다.

개인적으로 아쉬운 점: 후반부에 루즈해지는 것은 여전한 것 같다. bjwj5505가 G를 짤 때 C나 F를 제대로 봤어야 했다.

개인적으로 잘한 점: 푼 문제는 둘 다 원샷으로 끝냈고, K를 빠르게 풀어낸 것도 잘한 것 같다. 뇌절은 안 한듯.


풀이는 추후에 추가.

'PS > 대회 후기' 카테고리의 다른 글

SCPC 2020 1차 예선 풀이  (0) 2020.08.22
APIO 2020 Open Contest  (0) 2020.08.22
나는 코더다 2018 송년대회  (5) 2018.12.22
NYPC 2018 본선 복기 및 후기  (8) 2018.10.28
구사과와 함께하는 인터넷 예선 연습 대회  (0) 2018.10.08

8문제. 대략 난이도순.


BOJ 17257 불확정성이 넘쳐흘러

BOJ 16998 It's a Mod, Mod, Mod World

BOJ 11691 LCM(i, j)

BOJ 4862 마지막 자리

BOJ 8293 Prime Prime Power

BOJ 14861 LCM 더하기

BOJ 16941번 꼬리별

BOJ 13714번 약수의 개수


'PS > 수학 계열 PS' 카테고리의 다른 글

BOJ #9267 - A + B (SEERC 2013)  (2) 2019.11.24
Project Euler 300문제  (0) 2019.11.23
7-8월의 수학 PS 일지  (3) 2019.08.28
Project Euler #645  (0) 2019.03.01
BOJ #16164 Möbius Madness  (0) 2019.01.07

10문제, 어느 정도 난이도 순.


BOJ 12972 GCD 테이블



BOJ 13890 Big Bang



BOJ 13891 Find C



BOJ 2814 최소인수



BOJ 14858 GCD 테이블과 연속 부분 수열



BOJ 10916 Xtreme GCD Sum



BOJ 9647 Exponential Towers



BOJ 7936 N의 존재



BOJ 17372 피보나치 수의 최대공약수의 합



BOJ 17405 피보나치 수의 최대공약수의 합처럼 보이지만...





'PS > 수학 계열 PS' 카테고리의 다른 글

Project Euler 300문제  (0) 2019.11.23
9월초 수학 PS 일지  (3) 2019.09.11
Project Euler #645  (0) 2019.03.01
BOJ #16164 Möbius Madness  (0) 2019.01.07
Tonelli-Shanks Algorithm  (1) 2019.01.01

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


문제


주어지는 양의 정수 \(N, L, K\)에 대하여, 다음 값을 \(10^9 + 7\)로 나눈 나머지를 구하시오. $$ \sum_{d=1}^N \mu (L \cdot d) \left \lfloor \frac{N}{d} \right \rfloor^K $$ 단, \(1 \le N \le 10^9\), \(1 \le L \le 10^{15}\), \(1 \le K \le 10^9\)이며, 시간 제한은 2.5초이다. 


스포 방지선

________________________________________________________________________________________________________________________

















________________________________________________________________________________________________________________________


풀이


우선 여기에 설명되어 있는 테크닉인 xudyh's sieve에 대한 이해가 필요하다. 

두 multiplicative function \(f(x)\)와 \(g(x)\)가 다음을 만족한다고 하자.

  • \(g(x)\)의 partial sum, 즉 \(\sum_{i=1}^n g(i)\)를 빠른 시간에 구할 수 있다. 
  • \((f * g)(x)\)의 partial sum, 즉 \(\sum_{i=1}^n (f*g)(i)\)를 빠른 시간에 구할 수 있다.
  • 단, 여기서 \((f*g)(x)\)는 \(f(x)\)와 \(g(x)\)의 Dirichlet Convolution이다. 
  • 여기서 빠른 시간은 상수 시간이나 로그 시간, 로그 제곱 시간 정도를 말한다.

그러면, 각 \(1 \le k \le n\)에 대해 \(\sum_{i=1}^{\lfloor n/k \rfloor} f(i)\)를 \(\mathcal{O}(n^{2/3})\) 정도 시간에 구할 수 있다. 

여기서 가능한 \(\lfloor \frac{n}{k} \rfloor\)의 값의 종류가 \(\mathcal{O}(\sqrt n)\)개임을 참고하자. 자세한 방법론과 구현 방법은 링크 참조. 


문제를 풀기 위해서, 먼저 \( \left \lfloor \frac{N}{d} \right \rfloor ^K\)라는 항부터 보자. \( \left \lfloor \frac{N}{d} \right \rfloor \)로 가능한 값의 종류가 \(\mathcal{O}(\sqrt N)\)개이므로, 그 값에 따라서 \(d\)를 \(\mathcal{O}(\sqrt N)\)개의 구간으로 묶어줄 수 있다. 예를 들어, \(lef_v \le d \le rig_v \iff \left \lfloor \frac{N}{d} \right \rfloor = v\)라고 하면, \(lef_v \le d \le rig_v\)에 대하여 \( \left \lfloor \frac{N}{d} \right \rfloor^K\)라는 값은 \(v^K\)로 일정하므로, \( \sum_{d=lef_v}^{rig_v} \mu(L \cdot d)\)를 구해주면 된다. 

또한, 모든 \(\left \lfloor \frac{N}{d} \right \rfloor\) 형태로 쓸 수 있는 자연수 \(v\)에 대하여, \(rig_v = \left \lfloor \frac{N}{v} \right \rfloor \)임을 쉽게 확인할 수 있다. 

그러므로 모든 \(\left \lfloor \frac{N}{d} \right \rfloor\) 형태로 쓸 수 있는 자연수 \(v\)에 대하여, \(rig_v\)와 \(lef_v-1\)은 모두  \(\left \lfloor \frac{N}{d} \right \rfloor\) 형태로 쓸 수 있다.


결국 남은 것은 각 \(1 \le k \le n\)에 대해 \(\sum_{d=1}^{\lfloor N/k \rfloor} \mu(L \cdot d) \)를 빠르게 구하는 것이다. 


우선 \(L\)이 square-free하지 않다면, \(\mu(L\cdot d)\)는 \(d\)와 관계 없이 무조건 \(0\)일 것이다. 

그렇지 않다면, \(\mu_L(d) = \frac{\mu (L \cdot d)}{\mu (L)}\)이라 하자. 이제 목표를 다시 쓰자. 

우리의 목표는 각 \(1 \le k \le n\)에 대해 \(\sum_{d=1}^{\lfloor N/k \rfloor} \mu_L(d) \)를 빠르게 구하는 것이다. 

 

Claim: \(\mu_L(x)\)는 multiplicative function이다. 

Proof of Claim: \(\text{gcd}(d,L) \neq 1\)이라면 \(\mu_L(d)=0\)이고, \(\text{gcd}(d,L)=1\)이면 \(\mu_L(d) = \mu(d)\)이다. 


이제 xudyh's sieve를 쓸 각을 재보자. \(g(x)=1\)을 선택한다. 이제 \( (\mu_L * g) (x)\)에 대해 생각해보자. 


\(p|L\)인 소수 \(p\)와 자연수 \(k\)에 대해, \((\mu_L * g)(p^k) = \sum_{i=0}^k \mu_L(p^i) = 1\)이다.

\(L\)의 약수가 아닌 소수 \(q\)와 자연수 \(k\)에 대해, \((\mu_L * g)(q^k) = \sum_{i=0}^k \mu_L(q^i) = 0\)이다.

\((\mu_L * g)(x)\)도 multiplicative function이니, 우리는 다음과 같은 결론을 얻을 수 있다. 

모든 소수 \(p\)에 대해 \(p|x \implies p|L\)이 성립하면, \((\mu_L * g)(x) = 1\)이고, 그렇지 않다면 \((\mu_L * g)(x)=0\)이다. 


결국 \(\sum_{i=1}^n (\mu_L * g)(i)\)는 \([1, n]\)에 속한 자연수 중 \(L\)의 소인수들로만 구성된 자연수의 개수이다. 

\(L\)의 소인수들로만 구성된 자연수는 생각보다 많지 않다. \(L \le 10^{15}\)이므로, \(L\)이 가질 수 있는 소인수의 개수도 적다. \([1,N]\)에 속한 자연수 중 \(L\)의 소인수들로만 구성된 자연수는 priority_queue 등을 통하여 모두 구할 수 있다. 

이제 \(\sum_{i=1}^n (\mu_L * g)(i)\)는 lower_bound 등의 함수를 통해 빠르게 구할 수 있다. 


결국 xudyh's sieve를 쓰기 위한 필요 조건이 모두 갖추어졌고, 구현만 하면 된다. 끝.



'PS > 수학 계열 PS' 카테고리의 다른 글

Project Euler 300문제  (0) 2019.11.23
9월초 수학 PS 일지  (3) 2019.09.11
7-8월의 수학 PS 일지  (3) 2019.08.28
Project Euler #645  (0) 2019.03.01
Tonelli-Shanks Algorithm  (1) 2019.01.01

Aim: Solve \(x^2 \equiv n \pmod p \) in \( \mathcal{O}( \log^2 p ) \).


Tonelli-Shanks Algorithm은 \(\mathbb{Z}/p\mathbb{Z}\)에서 이차방정식을 빠르게 풀어내는 알고리즘이다.


Step 1. \(n\)이 \(\text{mod   } p\)에서 quadratic residue인지 확인한다. \( n^{\frac{p-1}{2}} \equiv 1 \pmod p\)인지 여부를 확인하면 된다.

Step 2. \(p-1 = Q \cdot 2^S\)인 홀수 \(Q\)와 음이 아닌 정수 \(S\)를 구한다. 

Step 3. \(\text{mod   } p\)에서 quadratic non-residue인 자연수 \(z\)를 하나 찾는다. 

Step 4. \(M=S\), \(c=z^Q\), \(t=n^Q\), \(R = n^{\frac{Q+1}{2}} \)라 하자. 전부 \(\text{mod   }p\)로 본다.

Step 5. 이제 다음을 반복한다.

        • \(t=0\)이면 \(0\)을 출력한다. \(t=1\)이면 \(R\)을 출력한다. 
        • 이제 \(0<i<M\)이면서 \(t^{2^i} \equiv 1 \pmod{p}\)를 만족하는 최소의 자연수 \(i\)를 찾는다.
        • \(b=c^{2^{M-i-1}}\)이라 두고, \(M \leftarrow i\), \(c \leftarrow b^2\), \(t \leftarrow tb^2\), \(R \leftarrow Rb\)로 값을 바꿔준다.

알고리즘을 통해 출력된 값이 \(r\)이라면, \(r\)와 \(-r\)이 \(x^2 \equiv n \pmod p\)의 두 근이 된다.

이 알고리즘을 활용하여, \(ax^2+bx+c \equiv n \pmod p\)도 쉽게 풀 수 있다.  


Claim 1: 이 알고리즘으로 출력된 값이 \(r\)이라면, \(r^2 \equiv n \pmod p\)가 성립한다. 

Proof of Claim 1: 우선 \(c^{2^{M-1}} = -1\), \(t^{2^{M-1}} = 1\), \(R^2 = tn\)이 항상 성립함을 보이자. 

초기에 \(c^{2^{M-1}} = z^{Q \cdot 2^{S-1}} = z^{\frac{p-1}{2}} = -1\)이다. (\(z\)는 quadratic non-residue니까)

초기에 \(t^{2^{M-1}} = n^{Q \cdot 2^{S-1}} = n^{\frac{p-1}{2}} = 1\)이다. (\(n\)은 quadratic residue니까)

초기에 \(R^2 = n^{Q+1} = n \cdot n^Q = nt\)가 성립한다. 초기에서는 세 식이 모두 성립함을 확인했다.


Step 5에서 \(M, c, t, R\)이 \(M', c', t', R'\)으로 바뀐다고 하자. 

\(c'^{2^{M'-1}} = b^{2^i} = c^{2^{M-1}} = -1\)이므로, 바뀐 뒤에도 첫번째 식이 성립한다. 

\(t'^{2^{M'-1}} = t^{2^{i-1}} \cdot b^{2^i} = -1 \cdot -1 = 1\)이므로, 바뀐 뒤에도 두번째 식이 성립한다. (why?)

\(R'^2 = R^2b^2 = tnb^2 = t'n\)이므로, 바뀐 뒤에도 마지막 식이 성립한다.  


\(t^{2^{M-1}} = 1\)이므로, Step 5에서 항상 조건을 만족하는 \(i\)를 찾을 수 있다.

또한, Step 5를 반복하면서 항상 \(M\)의 값이 감소하게 된다. \(M'=i<M\)이기 때문이다. 

그러다가 \(M=1\)이 되면 그 순간 \(t=1\)임이 보장되며, \(R^2 = tn = n\)이므로 \(R\)은 합동방정식의 해가 된다. 


Claim 2: 이 알고리즘은 \(\mathcal{O}(\log^2 p) \)번의 사칙연산을 필요로 한다. 

Proof of Claim 2: Step 1, 2, 4는 로그 시간에 할 수 있음이 자명하다. 

Step 3는 그냥 \(z\)를 랜덤하게 잡아보면서 quadratic non-residue가 나오면 terminate하는 방식으로 \(z\)를 찾아도 된다. 

\(z\)가 quadratic non-residue일 확률이 \(\frac{1}{2}\) 정도이므로, 평균 2번 정도의 시행을 필요로 한다. 

Step 5의 loop는 최대 \(S = \mathcal{O}(\log p)\)번 돌아가고, 각 loop에서 소모하는 연산 횟수는 \(\mathcal{O}(M)\)번, 즉 \(\mathcal{O}(\log p)\)번이다.

그러므로 Step 5에서는 최대 \(\mathcal{O}(\log^2 p)\)번의 사칙연산을 활용한다. 

시간복잡도를 \(\mathcal{O}(S^2)\)로도 쓸 수 있으며, \(S\)가 작은 경우도 많다는 점도 주목하자.


이제 Hensel's Lemma를 걸어주면, \(x^2 \equiv n \pmod{p^k}\)도 빠르게 풀어줄 수 있다. 


활용문제: Project Euler 216, 437



'PS > 수학 계열 PS' 카테고리의 다른 글

Project Euler 300문제  (0) 2019.11.23
9월초 수학 PS 일지  (3) 2019.09.11
7-8월의 수학 PS 일지  (3) 2019.08.28
Project Euler #645  (0) 2019.03.01
BOJ #16164 Möbius Madness  (0) 2019.01.07