QOJ.ac

QOJ

ID题目提交者结果用时内存语言文件大小提交时间测评时间
#35757#1860. Historic Breakthroughhos_lyricTL 22ms3824kbC++1410.9kb2022-06-18 18:45:122022-06-18 18:45:14

Judging History

你现在查看的是最新测评结果

  • [2023-08-10 23:21:45]
  • System Update: QOJ starts to keep a history of the judgings of all the submissions.
  • [2022-06-18 18:45:14]
  • 评测
  • 测评结果:TL
  • 用时:22ms
  • 内存:3824kb
  • [2022-06-18 18:45:12]
  • 提交

answer

#pragma GCC optimize ("Ofast")
#pragma GCC optimize ("unroll-loops")

#include <cassert>
#include <cmath>
#include <cstdint>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <algorithm>
#include <bitset>
#include <complex>
#include <deque>
#include <functional>
#include <iostream>
#include <map>
#include <numeric>
#include <queue>
#include <set>
#include <sstream>
#include <string>
#include <unordered_map>
#include <unordered_set>
#include <utility>
#include <vector>

using namespace std;

using Int = long long;

template <class T1, class T2> ostream &operator<<(ostream &os, const pair<T1, T2> &a) { return os << "(" << a.first << ", " << a.second << ")"; };
template <class T> void pv(T a, T b) { for (T i = a; i != b; ++i) cerr << *i << " "; cerr << endl; }
template <class T> bool chmin(T &t, const T &f) { if (t > f) { t = f; return true; } return false; }
template <class T> bool chmax(T &t, const T &f) { if (t < f) { t = f; return true; } return false; }


constexpr unsigned __int128 toUInt128(const char *s) {
  unsigned __int128 x = 0;
  for (; *s; ++s) x = x * 10 + (*s - '0');
  return x;
}
constexpr __int128 toInt128(const char *s) {
  if (*s == '-') return -toInt128(s + 1);
  __int128 x = 0;
  for (; *s; ++s) x = x * 10 + (*s - '0');
  return x;
}
unsigned __int128 inUInt128() {
  static char buf[40];
  scanf("%s", buf);
  return toUInt128(buf);
}
__int128 inInt128() {
  static char buf[40];
  scanf("%s", buf);
  return toInt128(buf);
}

void out(unsigned __int128 x) {
  static char buf[40];
  int len = 0;
  do { buf[len++] = '0' + static_cast<int>(x % 10); } while (x /= 10);
  for (int i = len; --i >= 0; ) putchar(buf[i]);
}
void out(__int128 x) {
  if (x < 0) {
    putchar('-');
    out(-static_cast<unsigned __int128>(x));
  } else {
    out(static_cast<unsigned __int128>(x));
  }
}
std::ostream &operator<<(std::ostream &os, unsigned __int128 x) {
  static char buf[40];
  int len = 0;
  do { buf[len++] = '0' + static_cast<int>(x % 10); } while (x /= 10);
  for (int i = len; --i >= 0; ) os << buf[i];
  return os;
}
std::ostream &operator<<(std::ostream &os, __int128 x) {
  if (x < 0) {
    os << '-' << -static_cast<unsigned __int128>(x);
  } else {
    os << static_cast<unsigned __int128>(x);
  }
  return os;
}


////////////////////////////////////////////////////////////////////////////////
// floor(a b / 2^128)
inline unsigned __int128 mul128High(unsigned __int128 a, unsigned __int128 b) {
  const unsigned __int128 a0 = static_cast<unsigned long long>(a), a1 = a >> 64;
  const unsigned __int128 b0 = static_cast<unsigned long long>(b), b1 = b >> 64;
  return ((((a0 * b0) >> 64) + static_cast<unsigned long long>(a0 * b1) + static_cast<unsigned long long>(a1 * b0)) >> 64) + ((a0 * b1) >> 64) + ((a1 * b0) >> 64) + a1 * b1;
}

// Hold y := (2^128 x) mod M
template <int ID> struct RMModInt128 {
  static unsigned __int128 M;
  static unsigned __int128 INV_M;  // M^-1 mod 2^128
  static unsigned __int128 TWO256;  // 2^256 mod M
  static void setM(unsigned __int128 m) {
    assert(m & 1); assert(1 <= m); assert(!(m >> 127));
    M = m;
    INV_M = (3 * M) ^ 2;
    for (int i = 0; i < 5; ++i) INV_M *= (2 - M * INV_M);
    TWO256 = (-M) % M;
    for (int i = 0; i < 128; ++i) TWO256 = ((TWO256 <<= 1) >= M) ? (TWO256 - M) : TWO256;
  }
  // (2^-128 a) mod M  (0 <= a < 2^128 m)
  static inline unsigned __int128 reduce(unsigned __int128 a) {
    const unsigned __int128 c = -mul128High(INV_M * a, M);
    return (c >= M) ? (c + M) : c;
  }
  // (2^-128 a b) mod M  (0 <= a b < 2^128 m)
  static inline unsigned __int128 mulReduce(unsigned __int128 a, unsigned __int128 b) {
    const unsigned __int128 c = mul128High(a, b) - mul128High(INV_M * (a * b), M);
    return (c >= M) ? (c + M) : c;
  }

  unsigned __int128 y;
  RMModInt128() : y(0U) {}
  RMModInt128(unsigned x_) : y(mulReduce(TWO256, x_ % M)) {}
  RMModInt128(unsigned long long x_) : y(mulReduce(TWO256, x_ % M)) {}
  RMModInt128(unsigned __int128 x_) : y(mulReduce(TWO256, x_ % M)) {}
  RMModInt128(int x_) : y(mulReduce(TWO256, ((x_ %= static_cast<__int128>(M)) < 0) ? (x_ + static_cast<__int128>(M)) : x_)) {}
  RMModInt128(long long x_) : y(mulReduce(TWO256, ((x_ %= static_cast<__int128>(M)) < 0) ? (x_ + static_cast<__int128>(M)) : x_)) {}
  RMModInt128(__int128 x_) : y(mulReduce(TWO256, ((x_ %= static_cast<__int128>(M)) < 0) ? (x_ + static_cast<__int128>(M)) : x_)) {}
  unsigned __int128 get() const { return reduce(y); }
  RMModInt128 &operator+=(const RMModInt128 &a) { y = ((y += a.y) >= M) ? (y - M) : y; return *this; }
  RMModInt128 &operator-=(const RMModInt128 &a) { y = ((y -= a.y) >= M) ? (y + M) : y; return *this; }
  RMModInt128 &operator*=(const RMModInt128 &a) { y = mulReduce(y, a.y); return *this; }
  RMModInt128 &operator/=(const RMModInt128 &a) { return (*this *= a.inv()); }
  template <class T> RMModInt128 pow(T e) const {
    if (e < 0) return inv().pow(-e);
    for (RMModInt128 a = *this, b = 1U; ; a *= a) { if (e & 1) { b *= a; } if (!(e >>= 1)) { return b; } }
  }
  RMModInt128 inv() const {
    unsigned __int128 a = M, b = reduce(reduce(y)); __int128 u = 0, v = 1;
    for (; b; ) { const unsigned __int128 q = a / b; const unsigned __int128 c = a - q * b; a = b; b = c; const __int128 w = u - static_cast<__int128>(q) * v; u = v; v = w; }
    assert(a == 1U); RMModInt128 r; r.y = u; r.y = (r.y >= M) ? (r.y + M) : r.y; return r;
  }
  RMModInt128 operator+() const { return *this; }
  RMModInt128 operator-() const { RMModInt128 a; a.y = y ? (M - y) : 0U; return a; }
  RMModInt128 operator+(const RMModInt128 &a) const { return (RMModInt128(*this) += a); }
  RMModInt128 operator-(const RMModInt128 &a) const { return (RMModInt128(*this) -= a); }
  RMModInt128 operator*(const RMModInt128 &a) const { return (RMModInt128(*this) *= a); }
  RMModInt128 operator/(const RMModInt128 &a) const { return (RMModInt128(*this) /= a); }
  template <class T> friend RMModInt128 operator+(T a, const RMModInt128 &b) { return (RMModInt128(a) += b); }
  template <class T> friend RMModInt128 operator-(T a, const RMModInt128 &b) { return (RMModInt128(a) -= b); }
  template <class T> friend RMModInt128 operator*(T a, const RMModInt128 &b) { return (RMModInt128(a) *= b); }
  template <class T> friend RMModInt128 operator/(T a, const RMModInt128 &b) { return (RMModInt128(a) /= b); }
  explicit operator bool() const { return y; }
  bool operator==(const RMModInt128 &a) const { return (y == a.y); }
  bool operator!=(const RMModInt128 &a) const { return (y != a.y); }
  friend std::ostream &operator<<(std::ostream &os, const RMModInt128 &a) { return os << get(); }
};
template <int ID> unsigned __int128 RMModInt128<ID>::M;
template <int ID> unsigned __int128 RMModInt128<ID>::INV_M;
template <int ID> unsigned __int128 RMModInt128<ID>::TWO256;
////////////////////////////////////////////////////////////////////////////////


using RMM128ForPrime = RMModInt128<-20220617>;

template <class T> vector<T> merge(const vector<T> &a, const vector<T> &b) {
  vector<T> c(a.size() + b.size());
  std::merge(a.begin(), a.end(), b.begin(), b.end(), c.begin());
  return c;
}

int bsf128(__int128 a) {
  const long long a64 = a;
  return a64 ? __builtin_ctzll(a64) : (64 + __builtin_ctzll(a >> 64));
}
__int128 gcd128(__int128 a, __int128 b) {
  if (a < 0) a = -a;
  if (b < 0) b = -b;
  if (a == 0) return b;
  if (b == 0) return a;
  const int s = bsf128(a | b);
  a >>= bsf128(a);
  do {
    b >>= bsf128(b);
    if (a > b) swap(a, b);
    b -= a;
  } while (b);
  return a << s;
}

// Checks if n is a prime using Miller-Rabin test
bool isPrime128(__int128 n) {
  if (n <= 1 || !(n & 1)) return (n == 2);
  RMM128ForPrime::setM(n);
  const int s = bsf128(n - 1);
  const __int128 d = (n - 1) >> s;
  // Based on conjectures in Zhang, Two kinds of strong pseudoprimes up to 10^36.
  for (const __int128 base : {2, 3, 5, 7, 11, 13, 17, 19, 23, 29,
                              31, 37, 41, 43, 47, 53, 59, 61, 67, 71}) {
    if (base >= n) continue;
    RMM128ForPrime a = RMM128ForPrime(base).pow(d);
    if (a.get() == 1 || static_cast<__int128>(a.get()) == n - 1) continue;
    bool ok = false;
    for (int i = 0; i < s - 1; ++i) {
      if (static_cast<__int128>((a *= a).get()) == n - 1) {
        ok = true;
        break;
      }
    }
    if (!ok) return false;
  }
  return true;
}

// Factorize n using Pollard's rho algorithm
vector<__int128> factorize128(__int128 n) {
  static constexpr int BLOCK = 256;
  if (n <= 1) return {};
  if (isPrime128(n)) return {n};
  if (!(n & 1)) {
    const int s = bsf128(n);
    return merge(vector<__int128>(s, 2), factorize128(n >> s));
  }
  RMM128ForPrime::setM(n);
  for (__int128 c0 = 2; ; ++c0) {
    const RMM128ForPrime c = c0;
    RMM128ForPrime x, y = 2, y0, z = 1;
    __int128 d = 1;
    for (int l = 1; d == 1; l <<= 1) {
      x = y;
      for (int i = 0; i < l; ++i) y = y * y + c;
      for (int i = 0; i < l; i += BLOCK) {
        y0 = y;
        for (int j = 0; j < BLOCK && j < l - i; ++j) {
          y = y * y + c;
          z *= (y - x);
        }
        if ((d = gcd128(z.y, n)) != 1) break;
      }
    }
    if (d == n) {
      for (y = y0; ; ) {
        y = y * y + c;
        if ((d = gcd128((y - x).y, n)) != 1) break;
      }
    }
    if (d != n) return merge(factorize128(d), factorize128(n / d));
  }
}

////////////////////////////////////////////////////////////////////////////////


__int128 solve(const __int128 M) {
// cerr<<"M = "<<M<<endl;
  vector<__int128> ps;
  
  // p (p-1)
  {
    __int128 lo = 1, hi = 2'000'000'000'000'000'000LL + 1;
    for (; lo + 1 < hi; ) {
      const __int128 mid = (lo + hi) / 2;
      ((mid * (mid - 1) <= M) ? lo : hi) = mid;
    }
    if (lo * (lo - 1) == M) {
// cerr<<"LINE "<<__LINE__<<" lo = "<<lo<<endl;
      if (isPrime128(lo)) {
        return lo;
      }
      ps = merge(factorize128(lo), factorize128(lo - 1));
      goto factored;
    }
  }
  
  ps = factorize128(M);
 factored:{}
// cerr<<"ps = ";pv(ps.begin(),ps.end());
  
  map<__int128, int> pes;
  for (const __int128 p : ps) {
    ++pes[p];
  }
  
  __int128 ans = 1;
  for (; !pes.empty(); ) {
    auto it = pes.end(); --it;
    const __int128 p = it->first;
    const int e = it->second;
    pes.erase(it);
    if (e > 0) {
      assert(e % 2 != 0);
      for (int i = 0; i < (e + 1) / 2; ++i) {
        ans *= p;
      }
      __int128 m = p - 1;
      for (auto &qf : pes) {
        for (; m % qf.first == 0; m /= qf.first) {
          --qf.second;
        }
        assert(qf.second >= 0);
      }
      assert(m == 1);
    }
  }
  return ans;
}

int main() {
  for (int numCases; ~scanf("%d", &numCases); ) for (; numCases--; ) {
    const __int128 M = 2 * inInt128();
    const __int128 ans = solve(M);
    out(ans);
    puts("");
  }
  return 0;
}

详细

Test #1:

score: 100
Accepted
time: 3ms
memory: 3808kb

input:

3
20
21
475750381222420656586462245096576000

output:

10
7
1497700821900508526

result:

ok 3 number(s): "10 7 1497700821900508526"

Test #2:

score: 0
Accepted
time: 22ms
memory: 3824kb

input:

51
348387408908517538156281238966340503
269891120302452381431351214335847781
747207543121624879797402427368860
500118165772005573992050274078796601
376563350255195175098956276198783051
855996192374691515214841787085600
491448606692239765059794615991064231
123619467864337410141102480048000000
7114827...

output:

834730386302688203
734698741393303847
38657773487574029
1000118158791255599
867828727636041299
41376351752391209
991411727479799147
819415677966571060
533472702079376326
419694411774324997
119851595982618319
24024477947405473
730267680763188643
269435681305057117
809387811759102827
29392724088327775...

result:

ok 51 numbers

Test #3:

score: -100
Time Limit Exceeded

input:

50
590955592280751522125185760551589472
450753984250583112023852704149662080
196704025354160782063198166237303808
382785853244719627595443384812477912
40522659517926585041466149305181616
26478235572251423131073298958930080
490320199321080674802144988571268192
110281951063110963040645709560838400
948...

output:


result: