저도 코드잼 D Large 해결... - 네타 주의 !

  • astein
    astein
    #include <iostream>
    #include <algorithm>
    #include <cstdio>
    #include <vector>
    
    using namespace std;
    
    // Maximum Number
    #define MAXNUM 1000000000
    // MAXNUM2 = MAXNUM/2
    #define MAXNUM2 500000000
    // SQRTMAXNUM = SQRT(MAXNUM)
    #define SQRTMAXNUM 31622
    // SQRTMAXNUM2 = SQRT(MAXNUM)/2
    #define SQRTMAXNUM2 15811
    // ARRAYSIZE = MAXNUM/2/32 + 1
    #define ARRAYSIZE MAXNUM/2/32 + 1
    #define MAXM 316227766
    
    #define getbit(n) (prime[((n)>>5)]&(1<<((n)&0x1F)))
    #define setbit(n) (prime[((n)>>5)]|=(1<<((n)&0x1F)))
    
    // isPrime() MACRO version
    #define isPrime(n) (( (((n)&1) && !(getbit((n)>>1))) || ((n)==2)) ? 1 : 0)
    
    unsigned int prime[ARRAYSIZE];
    int p[51000000];
    // MINimum MULtiplier
    int minmul[20001];
    int pn;
    
    // 소수 구하는 부분
    void initPrime() {
      setbit(0);
      for (int i = 1 ; i <= SQRTMAXNUM2 ; i++) {
        if (getbit(i) == 0) {
          for ( int j = ((i*i + i)<<1) ; j <= MAXNUM2 ; j += ((i<<1)|1))
          setbit(j);
        }
      }
      p[0] = 2;
      pn = 1;
      for (int i = 3; i < MAXNUM; i += 2) {
        if ( getbit((i - 1) / 2) == 0 ) {
          p[pn++] = i;
        }
      }
    }
    
    // 약수가 몇개인지 구하는 부분
    int numDivisor(long long n) {
      int ret = 1;
      int maxi = sqrt(n * 1.0);
      for (int i = 0; i <= maxi; ++i) {
        if (n % p[i] == 0) {
          int pp = 1;
          while (n % p[i] == 0) {
            n /= p[i];
            pp++;
          }
          ret *= pp;
        }
      }
      if (n != 1) ret *= 2;
      return ret;
    }
    
    // a 이상 b 이하의 소수의 개수를 구한다.
    // 이미 p 배열에는 10^9 이하의 모든 소수가 순서대로 저장되어 있기 때문에
    // upper_bound를 사용하면 ~O(lg n)만에 구할 수 있다.
    int numPrime(int a, int b) {
      return max(0, upper_bound(p, p + pn, b) - upper_bound(p, p + pn, a - 1));
    }
    
    // mpos는 앞으로 사용할 제일 작은 prime number의 index를 가지고 있다.
    // nd는 약수의 개수를 의미한다.
    long long find_sol(long long n, int mpos, int nd) {
      // n, mpos에 관계없이 nd = 1인 경우에는 답이 존재한다.
      if (nd == 1) return 1;
      // nd가 2인 경우에는 p[mpos] 이상 n 이하의 모든 소수가 답이 될 수 있다.
      if (nd == 2) return numPrime(p[mpos], (int)n);
      // p[mpos - 1] 이상의 소수를 사용하는 답은 존재하지 않는다.
      if (mpos == pn) return 0;
    
      long long ret = 0;
      // 만약 nd = 12인 경우 a^11, a^5 * b, a^3 * b^2, a^2 * b * c 등의 형태가 존재하는데,
      // a < b < c일 때 위의 4개의 경우 중 어떤 것을 고르더라도, a^4 보다는 크게 되기 때문에
      // x^4 <= n를 만족하는 x에 대해서 답을 구한다.
      int maxi = (int)(pow(n, 1.0 / minmul[nd]) + 1);
      for (int i = mpos; p[i] <= maxi; ++i) {
        long long t = 1;
        for (int j = 2; j <= nd; ++j) {
          // 주의: t * p[i]는 overflow가 날 수 있음
          if (t > n / p[i]) {
            break;
          }
          t = t * p[i];
          if (nd % j == 0) ret += find_sol(n / t, i + 1, nd / j);
        }
      }
      return ret;
    }
    
    int main() {
      initPrime();
      for (int i = 2; i <= 20000; ++i) {
        minmul[i] = i - 1;
        for (int j = 2; j * j <= i; ++j) {
          if (i % j == 0) {
            if (minmul[i] > minmul[j] + minmul[i / j]) {
              minmul[i] = minmul[j] + minmul[i / j];
            }
          }
        }
      }
    
      freopen("d-large-practice.in", "r", stdin);
      FILE *out = fopen("d-large-practice.out", "w");
      int T;
      cin >> T;
    
      long long n, m;
      for (int cn = 1; cn <= T; ++cn) {
        cin >> n >> m;
        int nd = numDivisor(n);
        // n = 1 또는 소수이거나, m^2이 n 이상이면 무조건 답은 0
        if (n == 1 || nd == 2 || m * m >= n) {
          printf("Case #%d: 0\n", cn);
          fprintf(out, "Case #%d: 0\n", cn);
        } else {
          // 입력 받은 m 이상의 제일 작은 소수의 위치를 찾음
          int mpos = lower_bound(p, p + pn, m) - p;
          long long ret = find_sol(n - 1, mpos, nd);
          printf("Case #%d: %I64d\n", cn, ret);
          fprintf(out, "Case #%d: %I64d\n", cn, ret);
        }
      }
    }
    

    저도 풀었습니다!

    기본 알고리즘은 아래 Andromeda님이 써 주신 내용과 비슷하고요.

    Large Data 1000개 전부 구하는데 30초정도 걸리네요.

    속도를 개선하는 요소는 크게 두 가지입니다.

    1) [A, B] 구간의 소수를 빨리 구한다.

    • 일단 N이 소수인 경우에는 항상 답이 0이니까, N이 두 개의 소수가 곱해진 형태의 합성수를 살펴봅시다. 이 경우 N이 10^9보다 작은 충분히 큰 경우의 수는 꽤 많이 존재합니다. 이를 일일이 세어본다면 8분의 시간으로는 어림도 없지요. 그래서 소수*소수 꼴에서 앞에 들어갈 수를 결정했다면 뒤에 들어갈 수의 범위가 나오게 됩니다. 여기에 들어갈 수 있는 모든 소수를 미리 계산해뒀다가, binary search를 응용해서 lg n 타임에 구해야 합니다. 참고로 대략 10^9 까지의 소수를 5천만개정도 있습니다.

    2) N 이하의 수 중에서 약수가 nd개인 경우가 있다.

    • 시간만 무한하다면야 모든 경우를 세는것이 제일 확실한 방법이지만, 우리에게는 8분의 여유밖에 없습니다. 따라서 답이 될 수 없는 경우는 탐색을 하지 않는 것이 필요하지요.
      D번 문제에서 큰 데이터를 테스트 하신 분들은 nd = 8 또는 16일 때 프로그램의 실행속도가 매우 느려지는 현상이 발생할 것입니다. N = 1억 정도이고 m = 2인 입력이 들어왔다면 소수 3개의 곱이 1억이 되지 않는 경우는 아주 많기 때문이지요. 위의 1)을 통해서 많이 개선이 되긴 하지만 여전히 느린것은 마찬가지입니다.
      답이 p1 * p2 * p3 (단 p1 < p2 < p3) 인 경우, p1에 2부터 차례대로 소수를 넣고 테스트하게 되는데, 여기에서 답이 존재하는 최대 p1이 몇인지 바로 판단할 수만 있다면 많은 시간을 아낄 수 있을 것입니다. p1 = 2나 3일 경우에는 답이 존재하겠지만 1000을 넘는 소수라면 p2, p3를 아무리 골라도 세 수의 곱은 1억을 넘게 될 테니까요.
      즉 p1 ^ 3 > N을 만족하는 p1은 그 다음수를 어떻게 고르더라도 해가 될 수 없다는 것을 알 수 있습니다. (물론 p1 ^ 3 <= N 을 만족한다고 해서 p1에 항상 답을 가진다는 보장은 없습니다.)
      nd = 8일때는 세제곱으로 테스트로 충분한 이유는 a^7, a^3 * b, a * b * c의 세 가지 형태 중에서 모든 소수에 a를 대입했을 때 세제곱이 나오기 때문이지요. 비슷한 예로 nd = 12일 때에는 a^11, a^5 * b, a^3 * b^2, a^2 * b * c가 있고, 이는 네제곱으로 테스트를 해서 a의 상한을 정할 수 있습니다.
      소스코드에서 minmul 배열은 nd가 주어졌을 때 몇제곱이 상한인지를 계산하고 있습니다. 그리고 실제로 k제곱 하게 되면 수행시간이 오래 걸리게 되므로 x ^ k <= n을 비교하는 대신 x <= n ^ (1 / k)를 비교합니다. 이렇게 하면 n과 k는 고정이기 때문에 매번 제곱연산을 할 필요가 없게 되지요.

      조만간 C번도 풀어보겠습니다.


    12년 전
6개의 댓글이 있습니다.
  • VOCList
    VOCList

    1빠


    12년 전 link
  • Andromeda
    Andromeda

    제 코드보다 3배쯤 빠르시네요. ㅎㄷㄷ
    많은 도움 되었습니다.


    12년 전 link
  • TurtleShip
    TurtleShip

    좋은 설명 감사드립니다. 혹시 소수 구하는 부분에서 setbit과 getbit에 대해서 부연 설명 가능하신가요?


    12년 전 link
  • TurtleShip
    TurtleShip

    astein님, astein님의 코드를 썼을 경우
    N = 6140400029153458
    M = 6140401
    에서 Floating point exception: 8 이라는 에러가 나요.

    디버깅 해보니, numDivisor(n)에서,

    maxi = 78360704 인데, 이 값은
    소수를 넣어놓은 배열 p의 크기 51000000 보다 커서 메모리 참조값이 어긋나서 에러가 나오는 거 같아요. 2가지를 알고 싶어요.
    1. 왜 배열 p의 크기를 51000000 로 하셨나요?
    2. large data 테스트 할 때 위의 에러가 발생 하셨을 것 같은데 어떻게 해결하셨나요?


    12년 전 link
  • Andromeda
    Andromeda
    1. 10^9보다 작은 소수가 50847534개 있어서 그래요.
    2. 추축해보면, "i <= maxi"를 "p[i] <= maxi"로 바꾸면 되겠네요. (실수로 코드를 구버전으로 올리셨나봐요.)

    12년 전 link
  • astein
    astein

    아... p[i] <= maxi... ㅠ


    12년 전 link
  • 정회원 권한이 있어야 커멘트를 다실 수 있습니다. 정회원이 되시려면 온라인 저지에서 5문제 이상을 푸시고, 가입 후 7일 이상이 지나셔야 합니다. 현재 문제를 푸셨습니다.