자기곱 (COI 2008)

백준 링크

---

This section is intentionally left blank.

---

굉장히 좋은 문제라고 생각한다. 처음에 필요한 관찰들도 좋았고, 구현도 잘만 하면 깔끔한 문제.

정수 \(x\)의 자릿수 곱을 \(p(x)\)라고 하자. 첫 번째 관찰은 \(p(x) < x\)라는 것.

따라서 자기곱 \(s(x) = xp(x)\)에서, \(p(x)^2 < xp(x) = s(x) < 10^{18}\)이므로 \(p(x) < 10^9\)이다.

두 번째 관찰은 \(p(x) < 10^9\) 범위에서 가능한 \(p(x)\)값의 개수가 굉장히 적다는 것이다. 자릿수는 1~9범위에 있기 때문에 

(\(A > 0\)이므로 0이 있는 경우는 제외한다) 가능한 소인수는 \(2,3,5,7\) 뿐이고, \(2^{e_1}3^{e_2}5^{e_3}7^{e_4} < X\)를 만족하는 \((e_1,e_2,e_3,e_4)\)의 개수를 \(K(X)\)라고 하면, 잘해야  \(\log_{2}X\cdot\log_{3}X\cdot\log_{5}X\cdot\log_{7}X = O(\lg^{4}X)\)개고, 실제로 백트래킹을 해보면 \(5194\)개라고 한다. \((B = 10^9)\)

이제 가능한 \(5194\)개의 \(p\)값을 모두 돌면서 구간 \(\left[\lceil \frac{A}{p} \rceil, \lfloor \frac{B}{p} \rfloor \right]\)에 있는 수 \(x\) 중 \(p(x) = p\)를 만족하는 \(x\)의 개수를 세자. 재귀함수를 이렇게 설계한다.

f(num_digits, prefix, low, high, target)

prefix : 지금까지 고정된 수의 prefix이다. 현재 탐색하는 범위는 [prefix000...0, prefix999...9]이다.

num_digits : leading zero를 포함하여 지금까지 탐색한 자리수이다. num_digits = 18일 때 종료.

low, high : 각각 \(\lceil \frac{A}{p} \rceil, \lfloor \frac{B}{p} \rfloor\)이다. 재귀에서 넘겨줄 때 변하지 않는다.

target : 현재 만들고자 하는 product

실제 구현에서는 prefix000...0과 prefix999...9를 쉽게 구하기 위해 pow10 = 1000...0이라는 인자를 하나 더 넘겨주는데, 취향껏.

만약 [prefix000...0, prefix999...9]가 [low, high]와 완전히 겹친다면 그 구간을 볼 필요가 없다.

[prefix000...0, prefix999...9]가 [low, high]에 완전히 포함되는 경우, 더 이상 prefix를 볼 필요가 없이 그냥 [0,999...9] 중 자릿수 곱이 target이 되는 수의 개수를 찾으면 된다. 이는 memoization을 이용해 O(1)에 처리 가능하다.

결국 돌아야 하는 재귀의 수는 \(\log_{10}B\)에 비례하게 된다. 시간 복잡도는 \(O(K(\sqrt{B})\log_{10}B)\)이다. \(K(X) = O(lg^{4}X)\)이므로, 좀 더 elementary한 함수로 표현하면 \(O(\lg^{5} B)\)이다. 물론 실제 값은 이보다 매우 작다.

#include <bits/stdc++.h> using namespace std; typedef long long ll; const int LG = 18; enum MAXP{ P2 = 30, P3 = 19, P5 = 13, P7 = 11}; ll T[ LG ][ P2 ][ P3 ][ P5 ][ P7 ]; ll garbage; int p[4] = {2,3,5,7}; int ft[10][4] = { {0,0,0,0},{0,0,0,0},{1,0,0,0},{0,1,0,0},{2,0,0,0}, {0,0,1,0},{1,1,0,0},{0,0,0,1},{3,0,0,0},{0,2,0,0} }; int e[4]; //range_count : [lo,hi]의 정수 중 digit product가 2^e[0] 3^e[1] 5^e[2] 7^e[3]과 같은 수의 개수 ll range_count(int dig, ll pfx, ll p10, ll lo, ll hi){ ll mfx = pfx + p10 - 1; if(pfx > hi || mfx < lo) return 0; if(dig == LG){ //base case return !(e[0] | e[1] | e[2] | e[3]); } int reusable = (pfx >= lo && mfx <= hi); //reusable ll &rf = (reusable ? T[dig][e[0]][e[1]][e[2]][e[3]] : garbage); if(reusable && rf >= 0) return rf; p10 /= 10; //decrease range ll ret = 0; for(int apnd = !!pfx; apnd < 10; apnd++){ bool ok = 1; for(int i = 0; i < 4; i++) ok &= ft[apnd][i] <= e[i]; if(!ok) continue; for(int i = 0; i < 4; i++) e[i] -= ft[apnd][i]; ret += range_count(dig+1, pfx + p10 * apnd, p10, lo, hi); for(int i = 0; i < 4; i++) e[i] += ft[apnd][i]; } if(reusable) rf = ret; return ret; } ll A, B; ll self_product; inline ll ceil(ll a, ll b){return (a+b-1)/b;} //dfs : exhaustive search for all possible products. const ll Ten18 = 1'000'000'000'000'000'000LL; void dfs(ll hi, ll now, int w){ //w : # of distinct primes if( now > 1e9 + 1 || now * now > hi ) return; //p(x)^2 < p(x) * x < hi if( w == 4 ){ //rec(digit number = 0, prefix = 0(empty), pow10, low, high) self_product += range_count(0, 0, Ten18, ceil(A,now), B/now); return; } //dfs1 : stop gaining current prime exp dfs(hi, now, w+1); //dfs2 : increase current prime exp by 1 ++e[w]; dfs(hi, now * p[w], w); --e[w]; } int main(){ cin >> A >> B; memset(T, -1, sizeof(T)); dfs(B,1,0); cout << self_product << '\n'; }