結果

問題 No.577 Prime Powerful Numbers
ユーザー qwewe
提出日時 2025-05-14 12:59:40
言語 C++17
(gcc 13.3.0 + boost 1.87.0)
結果
AC  
実行時間 1,118 ms / 2,000 ms
コード長 12,034 bytes
コンパイル時間 1,807 ms
コンパイル使用メモリ 99,244 KB
実行使用メモリ 7,844 KB
最終ジャッジ日時 2025-05-14 13:00:47
合計ジャッジ時間 6,181 ms
ジャッジサーバーID
(参考情報)
judge5 / judge2
このコードへのチャレンジ
(要ログイン)
ファイルパターン 結果
sample AC * 1
other AC * 10
権限があれば一括ダウンロードができます
コンパイルメッセージ
main.cpp: In function ‘long long int integer_kth_root(long long int, int)’:
main.cpp:123:45: warning: integer overflow in expression of type ‘u128’ {aka ‘__int128’} results in ‘0x7fffffffffffffffffffffffffffffff’ [-Woverflow]
  123 |             u128 MAX_VAL = (u128(1) << 127) - 1;
      |                            ~~~~~~~~~~~~~~~~~^~~
main.cpp:152:46: warning: integer overflow in expression of type ‘u128’ {aka ‘__int128’} results in ‘0x7fffffffffffffffffffffffffffffff’ [-Woverflow]
  152 |              u128 MAX_VAL = (u128(1) << 127) - 1;
      |                             ~~~~~~~~~~~~~~~~~^~~
main.cpp: In function ‘int main()’:
main.cpp:258:50: warning: integer overflow in expression of type ‘u128’ {aka ‘__int128’} results in ‘0x7fffffffffffffffffffffffffffffff’ [-Woverflow]
  258 |                  u128 MAX_VAL = (u128(1) << 127) - 1; // Max value for signed __int128
      |                                 ~~~~~~~~~~~~~~~~~^~~
main.cpp:277:55: warning: integer overflow in expression of type ‘u128’ {aka ‘__int128’} results in ‘0x7fffffffffffffffffffffffffffffff’ [-Woverflow]
  277 |                       u128 MAX_VAL = (u128(1) << 127) - 1;
      |                                      ~~~~~~~~~~~~~~~~~^~~

ソースコード

diff #

#include <iostream>
#include <vector>
#include <cmath>
#include <numeric>
#include <algorithm> // For std::min, std::reverse
#include <string> // For optional printing of __int128

// Use __int128_t for calculations involving potentially large numbers up to 10^18
// Note: __int128_t is a GCC/Clang extension. Make sure the compiler supports it.
using u128 = __int128_t;

// Modular multiplication: (a * b) % mod using binary multiplication (Russian peasant multiplication)
// This prevents overflow when a * b might exceed maximum __int128 value, provided mod is within range.
u128 mod_mul(u128 a, u128 b, u128 mod) {
    u128 res = 0;
    a %= mod;
    while (b > 0) {
        if (b & 1) res = (res + a) % mod;
        // Double 'a' and take modulo. Check for potential overflow if needed, but typically 2*a is fine if a < mod.
        a = (a * 2) % mod; 
        b >>= 1; // Halve 'b'
    }
    return res;
}

// Modular exponentiation: (base^exp) % mod
// Uses modular multiplication to prevent intermediate overflows.
u128 mod_pow(u128 base, u128 exp, u128 mod) {
    u128 res = 1;
    base %= mod;
    while (exp > 0) {
        if (exp & 1) res = mod_mul(res, base, mod);
        base = mod_mul(base, base, mod);
        exp >>= 1;
    }
    return res;
}

// Miller-Rabin primality test
// Efficiently checks if a number `n` is likely prime.
// Deterministic for n up to certain bounds depending on the chosen bases.
bool is_prime(long long n) {
    if (n < 2) return false; // Primes are >= 2
    if (n == 2 || n == 3) return true; // 2 and 3 are prime
    if (n % 2 == 0 || n % 3 == 0) return false; // Check divisibility by 2 and 3

    // Write n-1 as d * 2^s
    u128 d = n - 1;
    int s = 0;
    while ((d & 1) == 0) {
        d >>= 1;
        s++;
    }

    // Use a set of bases known to be sufficient for deterministic Miller-Rabin test up to large numbers.
    // The set {2, 3, 5, 7, 11, 13, 17, 19, 23} is sufficient for n < 3,825,123,056,546,413,051 (~3.8 * 10^18),
    // which covers the input range N <= 10^18.
    int bases[] = {2, 3, 5, 7, 11, 13, 17, 19, 23}; 
    
    for (long long a_ll : bases) {
         u128 a = a_ll;
         if (n == a) return true; // If n itself is one of the small prime bases
         if (a >= n) continue; // Base must be < n
        
        u128 x = mod_pow(a, d, n); // Compute a^d % n
        if (x == 1 || x == n - 1) continue; // Passes the first check
        
        bool maybe_prime = false; // Flag to track if n passes subsequent checks for this base
        for (int r = 1; r < s; ++r) {
            x = mod_mul(x, x, n); // Compute x^2 % n
            if (x == n - 1) { // If x == n-1, n might be prime for this base
                maybe_prime = true;
                break;
            }
            // If x becomes 1, then n is composite (found non-trivial square root of 1 modulo n)
             if (x == 1) { 
                 return false;
             }
        }
        // If after all squaring steps, x is not n-1, then n is composite
        if (!maybe_prime) return false; 
    }
    return true; // If n passes checks for all bases, it's highly likely prime (deterministic for this range)
}

// Function to compute integer b-th root of Y
// Returns q such that q^b = Y, or 0 if no such integer q >= 2 exists.
long long integer_kth_root(long long Y, int b) {
    if (Y <= 1) return 0; // Prime powers must be > 1
    if (b == 1) return Y; // For b=1, root is Y itself. Primality check handled elsewhere.

    long long low = 2, high = Y;
    // Optimize search range for high using properties of powers:
    if (b >= 60) high = 3; // Since 2^60 > 10^18, root can only be 2. Upper bound 3 is safe.
    else if (b >= 40) high = std::min((long long)4, Y); // e.g. 4^40 > 10^24
    else if (b >= 30) high = std::min((long long)15, Y); // e.g. 15^30 > 10^35
    else if (b == 2) high = std::min((long long)1000000000 + 7, Y); // Upper bound approx sqrt(10^18)
    else if (b == 3) high = std::min((long long)1000000 + 7, Y); // Upper bound approx cbrt(10^18)
    else {
        // For other b, compute approximate root using floating point and add a small buffer.
        // Using long double for better precision.
         long double root_approx = pow((long double)Y, 1.0L / b);
         high = std::min((long long)(root_approx + 2.0L), Y); // Add buffer and ensure <= Y
    }
     
    if (low > high) return 0; // Search range is invalid

    long long ans = 0;

    // Binary search for the integer root
    while(low <= high) {
        long long mid = low + (high - low) / 2;
        if (mid <= 1) { // Ensure mid >= 2 for prime power base
             low = 2; continue; // Adjust lower bound and continue search
        }

        u128 val = 1;
        bool overflow = false;
        // Calculate mid^b safely using __int128_t
        for(int i = 0; i < b; ++i) {
            // Check for potential overflow before multiplication
            // Max value for signed __int128 is approximately (1 << 127) - 1
            u128 MAX_VAL = (u128(1) << 127) - 1; 
            // If val > MAX_VAL / mid, then val * mid would overflow
            if (val > MAX_VAL / mid) { 
                  overflow = true;
                  break;
            }
            val *= mid; // Perform multiplication
            // If current power value exceeds Y, it's too large
            if (val > Y) { 
                overflow = true; // Consider this as "overflowing" the target Y
                break;
            }
        }

        if (overflow) { // If mid^b > Y or overflowed __int128
            high = mid - 1; // Root must be smaller
        } else if (val == Y) {
            ans = mid; // Found exact integer root
            break; 
        } else { // val < Y
            low = mid + 1; // Root must be larger
        }
    }
    
    // Final check to ensure the found root `ans` indeed satisfies ans^b == Y
    if (ans >= 2) { 
        u128 check_val = 1;
        bool ok = true;
        for(int i=0; i<b; ++i) {
             u128 MAX_VAL = (u128(1) << 127) - 1;
             if (check_val > MAX_VAL / ans) { ok = false; break; } // Overflow check
            check_val *= ans;
            if (check_val > Y) { ok = false; break; } // Exceeds Y check
        }
        if (ok && check_val == Y) return ans; // Confirmed correct root
        else return 0; // Found value doesn't satisfy check, return 0
    }
    
    return 0; // No integer root >= 2 found
}


// Check if Y can be represented as q^b where q is a prime and b >= 1
bool is_prime_power(long long Y) {
    if (Y <= 1) return false; // Prime powers must be > 1
    
    // Check the b=1 case first: Is Y itself a prime?
    if (is_prime(Y)) {
        return true;
    }

    // Check for cases b >= 2
    int max_b = 0;
    // Calculate max possible exponent b such that 2^b <= Y
    if (Y >= 4) max_b = floor(log2((long double)Y)); // Need Y>=4 because 2^2=4 is smallest power with b>=2
    else return false; // If Y=2 or Y=3, they are primes, handled above. If Y=1, returned false.

    // Since N <= 10^18, Y < 10^18. log2(10^18) is approx 59.79. Max exponent is 60.
    if (max_b > 60) max_b = 60; 

    // Iterate through possible exponents b from 2 up to max_b
    for (int b = 2; b <= max_b; ++b) { 
        long long q = integer_kth_root(Y, b); // Find the integer b-th root
        // If integer root q >= 2 exists and q is prime
        if (q >= 2 && is_prime(q)) {
            // The integer_kth_root function ensures q^b == Y. If it returns q >= 2, it found a valid root.
             return true;
        }
    }

    return false; // No representation as q^b with b >= 1 found.
}

// Sieve primes up to a bound B. Using B=10^6 which is generally fast enough.
const int B = 1000000; 
std::vector<int> primes; // Stores primes up to B
bool is_composite[B + 1]; // Sieve array

// Sieve of Eratosthenes implementation
void sieve() {
    std::fill(is_composite, is_composite + B + 1, false);
    is_composite[0] = is_composite[1] = true; // 0 and 1 are not prime
    for (int p = 2; p * p <= B; ++p) {
        if (!is_composite[p]) { // If p is prime
            for (int i = p * p; i <= B; i += p) // Mark multiples of p as composite
                is_composite[i] = true;
        }
    }
    // Collect all primes found into the `primes` vector
    for (int p = 2; p <= B; ++p) {
        if (!is_composite[p]) {
            primes.push_back(p);
        }
    }
}


int main() {
    // Faster I/O
    std::ios_base::sync_with_stdio(false);
    std::cin.tie(NULL);

    // Precompute primes up to B using the sieve
    sieve();

    int Q; // Number of test cases
    std::cin >> Q;
    while (Q--) {
        long long N; // Input number for the current test case
        std::cin >> N;

        bool found = false; // Flag to indicate if a solution is found

        // Handle small N values. Prime powers are >= 2. Smallest sum is 2+2=4.
        if (N <= 3) {
            std::cout << "No\n";
            continue;
        }
        
        // Optimization based on N's parity:
        // If N is odd, one term must be even, the other odd.
        // The only even prime power is 2^a (a >= 1).
        // So N must be of the form 2^a + q^b where q is an odd prime.
        if (N % 2 != 0) { 
            u128 X = 2; // Start with 2^1
            for (int a = 1; ; ++a) {
                 if (X >= N) break; // If 2^a >= N, stop
                 long long Y = N - (long long)X; // Calculate the remaining part
                 // Check if Y is a prime power (q^b)
                 if (Y > 0 && is_prime_power(Y)) {
                    found = true; // Solution found
                    break;
                 }
                 
                 // Calculate next power of 2 (2^(a+1)) safely
                 u128 MAX_VAL = (u128(1) << 127) - 1; // Max value for signed __int128
                 if (X > MAX_VAL / 2) { break; } // Check for overflow before X*2
                 X *= 2;
            }
        } else { // N is even
             // N can be sum of two odd prime powers, or one/two even prime powers (powers of 2).
             // Iterate through primes p <= B (including p=2)
             for (int p : primes) {
                 u128 X = p; // Start with p^1
                 for (int a = 1; ; ++a) {
                      if (X >= N) break; // If p^a >= N, stop for this prime p
                      long long Y = N - (long long)X; // Calculate remaining part
                      // Check if Y is a prime power (q^b)
                      if (Y > 0 && is_prime_power(Y)) {
                         found = true; // Solution found
                         break;
                      }
                      
                      // Calculate next power p^(a+1) = X * p safely
                      u128 MAX_VAL = (u128(1) << 127) - 1; 
                      // Check for potential overflow before X * p
                      // Need p > 1 check because p=1 is not prime (though not in `primes` vector)
                      if (p > 1 && X > MAX_VAL / p) { 
                          break; // Overflow would occur
                      }
                      
                      X *= p; // Compute next power
                      // Safety check, should not happen for p > 0
                      if (X == 0 && p > 0) break; 
                 }
                 if (found) break; // If found for current p, move to next test case
             }
             // Note: This approach implicitly assumes that if a solution exists,
             // at least one of the primes p or q is <= B. This is a common simplification
             // in competitive programming unless explicitly stated otherwise or complex cases are needed.
             // Handling the case where both p, q > B would require additional logic, possibly becoming much harder.
        }

        // Output the result for the current test case
        if (found) {
            std::cout << "Yes\n";
        } else {
            std::cout << "No\n";
        }
    }
    return 0;
}
0