Felix's Library

This documentation is automatically generated by online-judge-tools/verification-helper

View the Project on GitHub fffelix-huang/CP-stuff

:heavy_check_mark: Integer Factorization (Pollard Rho 質因數分解)
(library/math/factorize.hpp)

質因數分解

long long a;
vector<long long> factors = factorize(a); // 回傳 a 的質因數 (排序)

因數

long long a;
vector<long long> divs = divisors(a); // 回傳 a 的所有因數 (排序)

題目

ABC293 F - Zero or One

Depends on

Verified with

Code

#pragma once
#include <vector>
#include <cassert>
#include <algorithm>
#include "../misc/type-traits.hpp"
#include "binary-gcd.hpp"
#include "safe-mod.hpp"
#include "is-prime.hpp"
#include "../random/rng.hpp"

namespace felix {

template<class T>
T pollard_rho(T n) {
	using U = internal::safely_multipliable_t<T>;
	if(n % 2 == 0) {
		return 2;
	}
	if(is_prime(n)) {
		return n;
	}
	while(true) {
		const T R = rng() % (n - 1) + 1;
		auto f = [&](T x) -> T {
			return internal::safe_mod<U>(U(x) * x + R, n);
		};
		T x = 1, y = 2, ys = 1, q = 1, g = 1;
		constexpr int m = 128;
		for(int r = 1; g == 1; r <<= 1) {
			x = y;
			for(int i = 0; i < r; i++) {
				y = f(y);
			}
			for(int k = 0; k < r && g == 1; k += m) {
				ys = y;
				for(int i = 0; i < std::min(m, r - k); i++) {
					y = f(y);
					q = internal::safe_mod<U>(U(q) * internal::safe_mod(x - y, n), n);
				}
				g = binary_gcd(q, n);
			}
		}
		if(g == n) {
			do {
				ys = f(ys);
				T x2 = internal::safe_mod(x - ys, n);
				g = binary_gcd(x2, n);
			} while(g == 1);
		}
		if(g != n) {
			return g;
		}
	}
	assert(false);
}

template<class T>
std::vector<T> factorize(T n) {
	if(n <= 1) {
		return {};
	}
	std::vector<T> res = {n};
	for(int i = 0; i < (int) res.size(); i++) {
		T p = pollard_rho(res[i]);
		if(p != res[i]) {
			res[i] /= p;
			res.push_back(p);
			i--;
		}
	}
	std::sort(res.begin(), res.end());
	return res;
}

template<class T>
std::vector<T> divisors(T n) {
	if(n == 0) {
		return {};
	}
	std::vector<std::pair<T, int>> v;
	for(auto p : factorize(n)) {
		if(v.empty() || v.back().first != p) {
			v.emplace_back(p, 1);
		} else {
			v.back().second++;
		}
	}
	std::vector<T> res;
	auto f = [&](auto f, int i, T x) -> void {
		if(i == (int) v.size()) {
			res.push_back(x);
			return;
		}
		for(int j = v[i].second; ; j--) {
			f(f, i + 1, x);
			if(j == 0) {
				break;
			}
			x *= v[i].first;
		}
	};
	f(f, 0, 1);
	std::sort(res.begin(), res.end());
	return res;
}

} // namespace felix
#line 2 "library/math/factorize.hpp"
#include <vector>
#include <cassert>
#include <algorithm>
#line 3 "library/misc/type-traits.hpp"
#include <numeric>

#include <type_traits>


namespace felix {

namespace internal {

#ifndef _MSC_VER
template<class T> using is_signed_int128 = typename std::conditional<std::is_same<T, __int128_t>::value || std::is_same<T, __int128>::value, std::true_type, std::false_type>::type;
template<class T> using is_unsigned_int128 = typename std::conditional<std::is_same<T, __uint128_t>::value || std::is_same<T, unsigned __int128>::value, std::true_type, std::false_type>::type;
template<class T> using make_unsigned_int128 = typename std::conditional<std::is_same<T, __int128_t>::value, __uint128_t, unsigned __int128>;
template<class T> using is_integral = typename std::conditional<std::is_integral<T>::value || is_signed_int128<T>::value || is_unsigned_int128<T>::value, std::true_type, std::false_type>::type;
template<class T> using is_signed_int = typename std::conditional<(is_integral<T>::value && std::is_signed<T>::value) || is_signed_int128<T>::value, std::true_type, std::false_type>::type;
template<class T> using is_unsigned_int = typename std::conditional<(is_integral<T>::value && std::is_unsigned<T>::value) || is_unsigned_int128<T>::value, std::true_type, std::false_type>::type;
template<class T> using to_unsigned = typename std::conditional<is_signed_int128<T>::value, make_unsigned_int128<T>, typename std::conditional<std::is_signed<T>::value, std::make_unsigned<T>, std::common_type<T>>::type>::type;
#else
template<class T> using is_integral = typename std::is_integral<T>;
template<class T> using is_signed_int = typename std::conditional<is_integral<T>::value && std::is_signed<T>::value, std::true_type, std::false_type>::type;
template<class T> using is_unsigned_int = typename std::conditional<is_integral<T>::value && std::is_unsigned<T>::value, std::true_type, std::false_type>::type;
template<class T> using to_unsigned = typename std::conditional<is_signed_int<T>::value, std::make_unsigned<T>, std::common_type<T>>::type;
#endif

template<class T> using is_signed_int_t = std::enable_if_t<is_signed_int<T>::value>;
template<class T> using is_unsigned_int_t = std::enable_if_t<is_unsigned_int<T>::value>;
template<class T> using to_unsigned_t = typename to_unsigned<T>::type;

template<class T> struct safely_multipliable {};
template<> struct safely_multipliable<short> { using type = int; };
template<> struct safely_multipliable<unsigned short> { using type = unsigned int; };
template<> struct safely_multipliable<int> { using type = long long; };
template<> struct safely_multipliable<unsigned int> { using type = unsigned long long; };
template<> struct safely_multipliable<long long> { using type = __int128; };
template<> struct safely_multipliable<unsigned long long> { using type = __uint128_t; };

template<class T> using safely_multipliable_t = typename safely_multipliable<T>::type;

}  // namespace internal


}  // namespace felix

#line 2 "library/math/binary-gcd.hpp"

namespace felix {

template<class T>
inline T binary_gcd(T a, T b) {
	if(a == 0 || b == 0) {
		return a | b;
	}
	int8_t n = __builtin_ctzll(a);
	int8_t m = __builtin_ctzll(b);
	a >>= n;
	b >>= m;
	while(a != b) {
		T d = a - b;
		int8_t s = __builtin_ctzll(d);
		bool f = a > b;
		b = f ? b : a;
		a = (f ? d : -d) >> s;
	}
	return a << (n < m ? n : m);
}

} // namespace felix

#line 2 "library/math/safe-mod.hpp"

namespace felix {

namespace internal {

template<class T>
constexpr T safe_mod(T x, T m) {
	x %= m;
	if(x < 0) {
		x += m;
	}
	return x;
}

} // namespace internal


} // namespace felix
#line 4 "library/math/pow-mod.hpp"

namespace felix {

namespace internal {

template<class T>
constexpr T pow_mod_constexpr(T x, long long n, T m) {
	using U = safely_multipliable_t<T>;
	if(m == 1) {
		return 0;
	}
	U r = 1, y = safe_mod(x, m);
	while(n) {
		if(n & 1) {
			r = (r * y) % m;
		}
		y = (y * y) % m;
		n >>= 1;
	}
	return r;
}

} // namespace internal


} // namespace felix

#line 4 "library/math/is-prime.hpp"

namespace felix {

namespace internal {

bool miller_rabin(long long n, std::vector<long long> x) {
	long long d = n - 1;
	d >>= __builtin_ctzll(d);
	for(auto a : x) {
		if(n <= a) {
			return true;
		}
		long long t = d;
		__uint128_t y = pow_mod_constexpr(a, d, n);
		while(t != n - 1 && y != 1 && y != n - 1ULL) {
			y = y * y % n;
			t <<= 1;
		}
		if(y != n - 1ULL && t % 2 == 0) {
			return false;
		}
	}
	return true;
}

} // namespace internal


bool is_prime(long long n) {
	if(n <= 1) {
		return false;
	}
	for(int p : {2, 3, 5, 7, 11, 13, 17, 19, 23, 29}) {
		if(n % p == 0) {
			return n == p;
		}
	}
	if(n < (1LL << 30)) {
		return internal::miller_rabin(n, {2, 7, 61});
	}
	return internal::miller_rabin(n, {2, 325, 9375, 28178, 450775, 9780504, 1795265022});
}

} // namespace felix

#line 2 "library/random/rng.hpp"
#include <chrono>

namespace felix {

inline unsigned long long rng() {
	static unsigned long long SEED = std::chrono::steady_clock::now().time_since_epoch().count();
	SEED ^= SEED << 7;
	SEED ^= SEED >> 9;
	return SEED;
}

} // namespace felix
#line 10 "library/math/factorize.hpp"

namespace felix {

template<class T>
T pollard_rho(T n) {
	using U = internal::safely_multipliable_t<T>;
	if(n % 2 == 0) {
		return 2;
	}
	if(is_prime(n)) {
		return n;
	}
	while(true) {
		const T R = rng() % (n - 1) + 1;
		auto f = [&](T x) -> T {
			return internal::safe_mod<U>(U(x) * x + R, n);
		};
		T x = 1, y = 2, ys = 1, q = 1, g = 1;
		constexpr int m = 128;
		for(int r = 1; g == 1; r <<= 1) {
			x = y;
			for(int i = 0; i < r; i++) {
				y = f(y);
			}
			for(int k = 0; k < r && g == 1; k += m) {
				ys = y;
				for(int i = 0; i < std::min(m, r - k); i++) {
					y = f(y);
					q = internal::safe_mod<U>(U(q) * internal::safe_mod(x - y, n), n);
				}
				g = binary_gcd(q, n);
			}
		}
		if(g == n) {
			do {
				ys = f(ys);
				T x2 = internal::safe_mod(x - ys, n);
				g = binary_gcd(x2, n);
			} while(g == 1);
		}
		if(g != n) {
			return g;
		}
	}
	assert(false);
}

template<class T>
std::vector<T> factorize(T n) {
	if(n <= 1) {
		return {};
	}
	std::vector<T> res = {n};
	for(int i = 0; i < (int) res.size(); i++) {
		T p = pollard_rho(res[i]);
		if(p != res[i]) {
			res[i] /= p;
			res.push_back(p);
			i--;
		}
	}
	std::sort(res.begin(), res.end());
	return res;
}

template<class T>
std::vector<T> divisors(T n) {
	if(n == 0) {
		return {};
	}
	std::vector<std::pair<T, int>> v;
	for(auto p : factorize(n)) {
		if(v.empty() || v.back().first != p) {
			v.emplace_back(p, 1);
		} else {
			v.back().second++;
		}
	}
	std::vector<T> res;
	auto f = [&](auto f, int i, T x) -> void {
		if(i == (int) v.size()) {
			res.push_back(x);
			return;
		}
		for(int j = v[i].second; ; j--) {
			f(f, i + 1, x);
			if(j == 0) {
				break;
			}
			x *= v[i].first;
		}
	};
	f(f, 0, 1);
	std::sort(res.begin(), res.end());
	return res;
}

} // namespace felix
Back to top page