题解 | 麦森数-NOIP2003普及组复赛D题

D-麦森数_NOIP2003普及组复赛

https://ac.nowcoder.com/acm/contest/231/D

题目描述

形如2^P-1的素数称为麦森数，这时P一定也是个素数。但反过来不一定，即如果P是个素数，2^P-1不一定也是素数。到1998年底，人们已找到了37个麦森数。最大的一个是 $P=3021377$ ，它有909526位。麦森数有许多重要应用，它与完全数密切相关。

任务：输入 $P（1000<P<3100000）$ ，计算2^P-1的位数和最后500位数字（用十进制高精度数表示）

输入描述:

只包含一个整数 $P（1000<P<3100000）$

输出描述:

第一行：十进制高精度数 $2^P-1$ 的位数。
第2-11行：十进制高精度数 $2^P-1$ 的最后500位数字。（每行输出50位，共输出10行，不足500位时高位补0）
不必验证 $2^P-1$ 与P是否为素数。

示例1

输入

1279

输出

386
00000000000000000000000000000000000000000000000000
00000000000000000000000000000000000000000000000000
00000000000000104079321946643990819252403273640855
38615262247266704805319112350403608059673360298012
23944173232418484242161395428100779138356624832346
49081399066056773207629241295093892203457731833496
61583550472959420547689811211693677147548478866962
50138443826029173234888531116082853841658502825560
46662248318909188018470682222031405210266984354887
32958028878050869736186900714720710555703168729087

解答

这道题的难点在于两个，一是求出2^N-1的位数，而是求出2^N-1的后500位。

实际上第二个比较容易实现，写个高精度外加快速幂就可以了。由于最后结果位数比较大，即使是效率很高的快速幂也会超时，所以我们要另想办法求2^N-1的位数。我们发现，2^N的一位不可能是0，2^N-1和2^N的位数实际上是一样的。而且，若某个数是10^K，那么这个数的位数是K+1。把2N写成(10^lg2)^N，求出 $lg2\times N+1$ 就是最后答案的位数。

 #include <cstdio>
 #include <cstring>
 #include <cmath>
 
 const int maxn = 1e6 + 5;
 
  int p, ans[maxn], x[maxn], al, xl, c[maxn], cl;
  
  inline void mul(int a[], int b[], int& al, int bl) {
      memset(c, 0, sizeof(c));
      cl = 1;
     for (int i = 1; i <= al; ++i)
         for (int j = 1; j <= bl; ++j) {
              if (i + j - 1 > 500) break;
             c[i + j - 1] += a[i] * b[j];
              if (c[i + j - 1] > 9) {
                 c[i + j] += c[i + j - 1] / 10;
                  c[i + j - 1] %= 10;
             }
          }
     cl = al + bl;
      while (!c[cl] && cl > 0) --cl;
     for (int i = 1; i <= cl; ++i) a[i] = c[i];
      al = cl;
 }
 
  inline void quickpow() {
      ans[1] = 1, al = 1;
      x[1] = 2, xl = 1;
      while (p) {
          if (p & 1) mul(ans, x, al, xl);
         mul(x, x, xl, xl);
          p >>= 1;
     }
  }
  
  int main() {
     scanf("%d", &p);
      printf("%d\n", (int)(log10(2) * p + 1));
      quickpow();
      --ans[1];
      for (int i = 1; i < al; ++i)
          if (ans[i] < 0) --ans[i + 1], ans[i] += 10;
    for (int i = 10; i >= 1; --i) {
         if (i != 10) putchar('\n');
         for (int j = 50; j >= 1; --j)
              printf("%d", ans[(i - 1) * 50 + j] ? ans[(i - 1) * 50 + j] : 0);
     }
      return 0;
  }

来源： Mr^Kevin