QOJ.ac

QOJ

IDProblemSubmitterResultTimeMemoryLanguageFile sizeSubmit timeJudge time
#328282#7881. Computational ComplexityqLCompile Error//C++2023.7kb2024-02-15 18:51:082024-02-15 18:51:09

Judging History

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

  • [2024-02-15 18:51:09]
  • 评测
  • [2024-02-15 18:51:08]
  • 提交

answer

#include <algorithm>
#include <bit>
#include <bitset>
#include <cassert>
#include <cmath>
#include <complex>
#include <cstdio>
#include <functional>
#include <locale>
#include <numeric>
#include <queue>
#include <string>
#include <utility>
#include <vector>

/**
 * Author: Cracker
 * Last update: 24.02.08
 */
namespace Default_Source {
#define UseStd(obj) using std::obj;
#define OPERATOR(T, U, OPT)                                       \
	friend constexpr T operator OPT(const T &lhs, const U &rhs) { \
		return T(lhs) OPT## = rhs;                                \
	}

namespace Type {
using i32 = int;
using u32 = unsigned;
using i64 = long long;
using u64 = unsigned long long;
#ifdef _GLIBCXX_USE_INT128
using i128 = __int128;
using u128 = unsigned __int128;
#endif
using f32 = float;
using f64 = double;
using f128 = long double;
template <typename T>
constexpr T Inf{std::numeric_limits<T>::max()};
template <typename T>
constexpr T HalfInf{std::numeric_limits<T>::max() / 2};
}  // namespace Type
using namespace Type;

namespace Utility {
using std::max;
using std::min;
using std::swap;
#ifndef likely
#define likely(x) __builtin_expect(!!(x), 1)
#endif
#ifndef unlikely
#define unlikely(x) __builtin_expect(!!(x), 0)
#endif
template <typename T>
constexpr void ckmin(T &lhs, const T &rhs) noexcept {
	lhs = lhs < rhs ? lhs : rhs;
}
template <typename T>
constexpr void ckmin(T &base, const std::initializer_list<T> &L) noexcept {
	for (const T &it : L) ckmin(base, it);
}
template <typename T>
constexpr void ckmax(T &lhs, const T &rhs) noexcept {
	lhs = lhs < rhs ? rhs : lhs;
}
template <typename T>
constexpr void ckmax(T &base, const std::initializer_list<T> &L) noexcept {
	for (const T &it : L) ckmax(base, it);
}
}  // namespace Utility
using namespace Utility;

namespace Traits {
template <typename _InputIterator>
using LegacyInputIterator = typename std::enable_if_t<std::is_convertible<
	typename std::iterator_traits<_InputIterator>::iterator_category,
	std::input_iterator_tag>::value>;
}

namespace Math {
UseStd(cos) UseStd(sin) UseStd(tan) UseStd(acos) UseStd(asin) UseStd(atan);
UseStd(abs) UseStd(ceil) UseStd(floor) UseStd(round);
template <typename T>
T ceil(T up, T down) {
	static_assert(std::is_arithmetic_v<T>);
	return (up + down - 1) / down;
}
template <typename T>
T floor(T up, T down) {
	static_assert(std::is_arithmetic_v<T>);
	return up / down;
}
template <typename T>
constexpr T Pi = acos<T>(-1);
}  // namespace Math

namespace IO {
namespace Buffer {
template <size_t BUFSIZE, typename ChatT = char>
class InBuffer {
#ifdef __unix__
#define fread fread_unlocked
#else
#endif
	FILE *fp;
	ChatT buf[BUFSIZE], *p, *ped;

   public:
	InBuffer(FILE *_fp = stdin) noexcept
		: fp{_fp}, buf{}, p{buf}, ped{buf + fread(buf, 1, BUFSIZE, fp)} {
		*ped = EOF;
	}
	ChatT getc() noexcept { return p == ped ? EOF : *p++; }
	ChatT pos() const noexcept { return *p; }
	void nxt() noexcept { p != ped && ++p; }
	bool ugetc(ChatT ch) noexcept { return p != buf && (*--p = ch, true); }
#undef fread
};
template <typename CharT>
class InBuffer<0, CharT> {
#ifdef __unix__
#define fgetc fgetc_unlocked
#else
#define fgetc _fgetc_nolock
#endif
	FILE *fp;

   public:
	InBuffer(FILE *_fp = stdin) noexcept : fp{_fp} {}
	CharT getc() const noexcept { return fgetc(fp); }
	CharT pos() const noexcept {
		CharT ch = fgetc(fp);
		return ungetc(ch, fp), ch;
	}
	void nxt() const noexcept { fgetc(fp); }
	bool ugetc(CharT ch) const noexcept { return ungetc(ch, fp); }
};
}  // namespace Buffer

static const std::locale Local{};
template <size_t BUFSIZE, typename CharT = char>
class InStream {
	Buffer::InBuffer<BUFSIZE, CharT> buf;

   public:
	InStream(FILE *_fp = stdin) noexcept : buf{_fp} {}
	template <typename _Rule>
	bool test(_Rule rule) const noexcept {
		return rule(buf.pos(), Local);
	}
	template <typename _Rule>
	void skip(_Rule rule) noexcept {
		while (rule(buf.pos(), Local)) buf.nxt();
	}
	CharT pos() noexcept { return buf.pos(); }
	void nxt() noexcept { buf.nxt(); }
	CharT getc() noexcept { return buf.getc(); }
	bool ugetc() noexcept { return buf.ugetc(); }
	bool eof() const noexcept { return buf.pos() == EOF; }
	operator bool() const noexcept { return buf.pos() != EOF; }

	template <typename T>
	std::enable_if_t<std::is_integral_v<T> && std::is_signed_v<T>, InStream &>
	operator>>(T &x) noexcept {
		x = 0, skip(std::isspace<CharT>);
		if (buf.pos() == '-')
			for (buf.nxt(); test(std::isdigit<CharT>);)
				x = x * 10 - (buf.getc() & 15);
		else
			for (; test(std::isdigit<CharT>);) x = x * 10 + (buf.getc() & 15);
		return *this;
	}
	template <typename T>
	std::enable_if_t<std::is_integral_v<T> && std::is_unsigned_v<T>, InStream &>
	operator>>(T &x) noexcept {
		x = 0, skip(std::isspace<CharT>);
		for (; test(std::isdigit<CharT>);) x = x * 10 + (buf.getc() & 15);
		return *this;
	}
	template <typename T>
	std::enable_if_t<std::is_floating_point_v<T>, InStream &> operator>>(
		T &x) noexcept {
		static T f;
		x = 0, skip(std::isspace<CharT>);
		if (buf.pos() == '-') {
			for (buf.nxt(); test(std::isdigit<CharT>);)
				x = x * 10 - (buf.getc() & 15);
			if (buf.pos() == '.')
				for (f = 0.1, buf.nxt(); test(std::isdigit<CharT>); f *= 0.1)
					x -= f * (buf.getc() & 15);
		} else {
			for (; test(std::isdigit<CharT>);) x = x * 10 + (buf.getc() & 15);
			if (buf.pos() == '.')
				for (f = 0.1, buf.nxt(); test(std::isdigit<CharT>); f *= 0.1)
					x += f * (buf.getc() & 15);
		}
		return *this;
	}
	InStream &operator>>(bool &boolen) noexcept {
		return skip(std::isspace<CharT>), boolen = buf.getc() == '1', *this;
	}
	InStream &operator>>(char &ch) noexcept {
		return skip(std::isspace<CharT>), ch = buf.getc(), *this;
	}
	InStream &operator>>(char *str) noexcept {
		for (skip(std::isspace<CharT>); test(std::isgraph<CharT>);)
			*str++ = buf.getc();
		return *str = '\0', *this;
	}
	InStream &operator>>(std::string &str) noexcept {
		str = "";
		for (skip(std::isspace<CharT>); test(std::isgraph<CharT>);)
			str += buf.getc();
		return *this;
	}

	InStream &operator>>(
		std::tuple<char *&, std::function<bool(char, std::locale)>>
			strtp) noexcept {
		auto &[str, rule] = strtp;
		for (skip(std::isspace<CharT>); test(rule);) *str++ = buf.getc();
		return *str = '\0', *this;
	}
	InStream &operator>>(
		std::tuple<std::string &, std::function<bool(char, std::locale)>>
			strtp) noexcept {
		auto &[str, rule] = strtp;
		str = "";
		for (skip(std::isspace<CharT>); test(rule);) str += buf.getc();
		return *this;
	}

	template <typename T1, typename T2>
	InStream &operator>>(std::pair<T1, T2> &ptt) noexcept {
		return *this >> ptt.first >> ptt.second;
	}
	template <typename... Args>
	InStream &operator>>(std::tuple<Args...> &tp) noexcept {
		return tp = reads<Args...>(), *this;
	}
	template <typename... Args>
	InStream &read(Args &...args) noexcept {
		return void(std::initializer_list<InStream>{*this >> args...}), *this;
	}
	template <typename T>
	T read() noexcept {
		static_assert(!std::is_pointer_v<T>);
		static T tmp;
		return read(tmp), tmp;
	}
	template <typename T>
	std::tuple<T> reads() noexcept {
		return std::make_tuple(read<T>());
	}
	template <typename T, typename... Args>
	typename std::enable_if<sizeof...(Args) != 0, std::tuple<T, Args...>>::type
	reads() noexcept {
		static T tmp;
		return read(tmp), std::tuple_cat(std::tuple<T>{tmp}, reads<Args...>());
	}
};
}  // namespace IO

namespace Arithmetic {
using std::complex;
template <typename T>
static constexpr T mul(T x, T y, T p) noexcept {
	if constexpr (sizeof(T) < sizeof(u64))
		return 1ll * x * y % p;
	else if constexpr (sizeof(T) < sizeof(u128))
		return i128(x) * y % p;
	else
		return x * y % p;
}
template <typename T, typename U>
static constexpr T qpow(T x, U pw) noexcept {
	T ret{1};
	for (; pw; pw >>= 1, x *= x)
		if (pw & 1) ret *= x;
	return ret;
}
template <typename T, T P>
class ModNum {
	static_assert(P >= 0);

   public:
	using un_type = std::make_unsigned_t<T>;
	using exp_type = ModNum<T, P - 1>;

   private:
	static_assert(std::is_integral_v<T>);
	un_type v;
	static un_type Mod;
	constexpr static un_type unorm(un_type x) noexcept {
		return x >= mod() ? x - mod() : x;
	}
	template <typename U>
	constexpr static un_type norm(U x) noexcept {
		return unorm(x < 0 ? x + mod() : x);
	}

   public:
	constexpr ModNum() noexcept : v{} {}
	template <typename U, typename = std::enable_if_t<std::is_integral_v<U>>>
	constexpr ModNum(const U &x) noexcept : v{norm(x % T(mod()))} {}
	constexpr ModNum(const ModNum &ot) noexcept = default;
	constexpr ModNum(ModNum &&ot) noexcept = default;
	constexpr ModNum &operator=(const ModNum &ot) noexcept {
		return v = ot.v, *this;
	}
	constexpr ModNum &operator=(ModNum &&ot) noexcept {
		return v = ot.v, *this;
	}
	constexpr explicit operator T() const noexcept { return v; }
	constexpr T val() const noexcept { return v; }
	constexpr static un_type mod() noexcept {
		if constexpr (P > 0)
			return P;
		else
			return Mod;
	}
	constexpr static void set(T _Mod) noexcept { Mod = _Mod; }
	constexpr ModNum &operator++() noexcept { return v = norm(++v), *this; }
	constexpr const ModNum operator++(i32) noexcept {
		un_type tmp = v;
		v = norm(++v);
		return ModNum(tmp);
	}
	constexpr ModNum operator-() const noexcept { return ModNum(-T(v)); }
	constexpr ModNum inv() const noexcept { return qpow(*this, mod() - 2); }
	constexpr ModNum &operator+=(ModNum ot) noexcept {
		return v = unorm(v + ot.v), *this;
	}
	constexpr ModNum &operator-=(ModNum ot) noexcept {
		return v = unorm(v + mod() - ot.v), *this;
	}
	constexpr ModNum &operator*=(ModNum ot) noexcept {
		return v = mul(v, ot.v, mod()), *this;
	}
	constexpr ModNum &operator/=(ModNum ot) noexcept {
		return *this *= ot.inv();
	}
	OPERATOR(ModNum, ModNum, +);
	OPERATOR(ModNum, ModNum, -);
	OPERATOR(ModNum, ModNum, *);
	OPERATOR(ModNum, ModNum, /);
	constexpr bool friend operator==(ModNum lhs, ModNum rhs) noexcept {
		return lhs.val() == rhs.val();
	}
	constexpr bool friend operator!=(ModNum lhs, ModNum rhs) noexcept {
		return lhs.val() != rhs.val();
	}
	template <typename U>
	constexpr ModNum pow(U pw) const noexcept {
		return qpow(*this, pw);
	}
	template <size_t BUFSIZE, typename CharT>
	constexpr friend IO::InStream<BUFSIZE, CharT> &operator>>(
		IO::InStream<BUFSIZE, CharT> &is, ModNum &x) noexcept {
		is.skip(std::isspace<CharT>);
		if (is.pos() == '-')
			for (is.nxt(); is.test(std::isdigit<CharT>);)
				x = x * 10 - (is.getc() & 15);
		else
			for (; is.test(std::isdigit<CharT>);) x = x * 10 + (is.getc() & 15);
		return is;
	}
};
template <typename T, T P>
ModNum<T, P>::un_type ModNum<T, P>::Mod;

template <i32 P>
using Z32 = ModNum<i32, P>;
template <i64 P>
using Z64 = ModNum<i64, P>;

/**
 * I don't know how it works.
 */
template <typename T>
constexpr T findPrimitiveRoot() {
	T i = 2;
	i32 k = std::countr_zero<typename std::make_unsigned_t<decltype(T::mod())>>(
		T::mod() - 1);
	for (;; ++i)
		if (qpow(i, (T::mod() - 1) / 2) != 1) break;
	return qpow(i, (T::mod() - 1) >> k);
}
template <typename T>
constexpr T primitiveRoot = findPrimitiveRoot<T>();
template <typename T>
constexpr T InvprimitiveRoot = qpow(primitiveRoot<T>, T::mod() - 2);

template <typename T>
std::vector<T> fac{1};
template <typename T>
std::vector<T> ifac{1};
template <typename T>
std::vector<T> iv{0};
template <typename T>
constexpr void pre_fac(u32 m) noexcept {
	if (fac<T>.size() < m) {
		u32 k = fac<T>.size();
		fac<T>.resize(m);
		for (; k < m; ++k) fac<T>[k] = fac<T>[k - 1] * k;
	}
}
template <typename T>
void pre_ifac(u32 m) noexcept {
	if (ifac<T>.size() < m) {
		u32 k = ifac<T>.size();
		pre_fac<T>(m), ifac<T>.resize(m);
		ifac<T>.back() = qpow(fac<T>.back(), T::mod() - 2);
		for (--m; m > k; --m) ifac<T>[m - 1] = ifac<T>[m] * m;
	}
}
template <typename T>
constexpr void pre_iv(u32 m) noexcept {
	if (iv<T>.size() < m) {
		u32 k = iv<T>.size();
		pre_ifac<T>(m), iv<T>.resize(m);
		for (; k < m; ++k) iv<T>[k] = ifac<T>[k] * fac<T>[k - 1];
	}
}
template <typename T>
constexpr T binom(u32 n, u32 m) noexcept {
	return fac<T>[n] * ifac<T>[m] * ifac<T>[n - m];
}
}  // namespace Arithmetic

namespace Poly {
enum Trans_mode { DIF = 0, DIT = 1 };
constexpr u32 toLim(u32 len) { return 1u << (32 - std::countl_zero(len)); }
using std::vector;

template <typename T>
class NTT {
	static vector<T> roots;
	static constexpr void InitRoot(u32 k) noexcept {
		u32 l = roots.size();
		roots.resize(k);
		for (u32 l2 = l << 1; l < k; l = l2, l2 <<= 1) {
			T p = Arithmetic::qpow(Arithmetic::primitiveRoot<T>,
								   (T::mod() - 1) / l2);
			for (u32 i = l; i < l2; ++ ++i)
				roots[i + 1] = (roots[i] = roots[i / 2]) * p;
		}
	}

   public:
	template <Trans_mode mode, typename Iterator>
	static constexpr void trans(Iterator p, u32 n) noexcept {
		InitRoot(n);
		if constexpr (mode == DIF)
			for (u32 l = n >> 1; l; l >>= 1)
				for (u32 t = 0; t < n; t += (l << 1))
					for (u32 k = 0; k < l; ++k) {
						T x = p[t + k], y = p[t + l + k];
						p[t + k] = x + y, p[t + l + k] = roots[k + l] * (x - y);
					}
		else {
			for (u32 l = 1; l < n; l <<= 1)
				for (u32 t = 0; t < n; t += (l << 1))
					for (u32 k = 0; k < l; ++k) {
						T x = p[t + k], y = roots[k + l] * p[t + l + k];
						p[t + k] = x + y, p[t + l + k] = x - y;
					}
			const T ivn = T(1) / n;
			for (u32 k = 0; k < n; ++k) p[k] *= ivn;
			std::reverse(p + 1, p + n);
		}
	}
};
template <typename T>
vector<T> NTT<T>::roots{1, 1};

template <typename Iterator>
constexpr void dot(Iterator f, Iterator g, u32 deg) {
	for (u32 k = 0; k < deg; ++k) f[k] *= g[k];
}
enum PolyMode { Coefficients = 0, PipValue = 1 };
template <typename T, typename _Trans = NTT<T>>
class poly : public vector<T> {
	static _Trans trans;
	PolyMode mode = PolyMode::Coefficients;

#define f (*this)
	static poly &times(poly &p, poly &q, u32 n) noexcept {
		return p.Dft(n), q.Dft(n), dot(p.begin(), q.begin(), n), p.Idft(n);
	}

   public:
	using vector<T>::vector;
	constexpr poly(u32 n, std::function<T(u32)> func) noexcept {
		redeg(n);
		for (u32 k = 0; k < n; ++k) f[k] = func(k);
	}
	constexpr void fill(std::function<T(u32)> func = [](u32 k) -> T {
		return {};
	}) noexcept {
		for (u32 k = 0; k < deg(); ++k) f[k] = func(k);
	}
	constexpr void fill(
		u32 n,
		std::function<T(u32)> func = [](u32 k) -> T { return {}; }) noexcept {
		redeg(n);
		for (u32 k = 0; k < deg(); ++k) f[k] = func(k);
	}
	constexpr u32 deg() const noexcept { return vector<T>::size(); }
	constexpr poly &redeg(u32 n) noexcept {
		return vector<T>::resize(n), *this;
	}
	constexpr poly slice(u32 length, u32 first = 0) const noexcept {
		return {vector<T>::begin() + first,
				vector<T>::begin() + Utility::min(deg(), first + length)};
	}
	constexpr poly &Dft(u32 n) noexcept {
		if (mode != PipValue)
			mode = PipValue, redeg(n),
			trans.template trans<DIF>(vector<T>::begin(), n);
		return *this;
	}
	constexpr poly &Idft(u32 n) noexcept {
		return mode = Coefficients,
			   trans.template trans<DIT>(vector<T>::begin(), n), *this;
	}
	constexpr poly operator-() const noexcept {
		poly ret = *this;
		for (T &it : ret) it = -it;
		return ret;
	}
	constexpr poly &operator+=(const poly &ot) noexcept {
		redeg(Utility::max(deg(), ot.deg()));
		for (u32 k = 0; k < ot.deg(); ++k) f[k] += ot[k];
		return *this;
	}
	constexpr poly &operator-=(const poly &ot) noexcept {
		redeg(Utility::max(deg(), ot.deg()));
		for (u32 k = 0; k < ot.deg(); ++k) f[k] -= ot[k];
		return *this;
	}
	constexpr poly &operator*=(const T &value) noexcept {
		for (T &it : *this) it *= value;
		return *this;
	}
	constexpr poly &operator*=(poly ot) noexcept {
		u32 degree = deg() + ot.deg() - 1;
		return times(*this, ot, toLim(degree)).redeg(degree);
	}
	OPERATOR(poly, poly, +);
	OPERATOR(poly, poly, -);
	OPERATOR(poly, poly, *);
	OPERATOR(poly, T, *);
	constexpr poly Dvt(u32 degree) const noexcept {
		poly ret(degree);
		for (u32 k = 1; k < Utility::min(deg(), degree + 1); ++k)
			ret[k - 1] = k * f[k];
		return ret;
	}
	constexpr poly Igt(u32 degree) const noexcept {
		poly ret(degree);
		Arithmetic::pre_iv<T>(degree);
		for (u32 k = Utility::min(deg(), degree - 1); k > 0; --k)
			ret[k] = Arithmetic::iv<T>[k] * f[k - 1];
		return ret;
	}
	template <typename U>
	constexpr poly &operator>>=(U t) noexcept {
		if (t >= deg()) return *this = {};
		for (u32 k = t; k < deg(); ++k) f[k - t] = f[k];
		return redeg(deg() - t);
	}
	template <typename U>
	constexpr poly &operator<<=(U t) noexcept {
		redeg(deg() + t);
		for (u32 k = t; k < deg(); ++k) f[k] = f[k - t];
		for (u32 k = 0; k < t; ++k) f[k] = 0;
		return redeg(deg() - t);
	}
	template <typename U>
	constexpr friend poly operator>>(poly lhs, U t) noexcept {
		return lhs >>= t;
	}
	template <typename U>
	constexpr friend poly operator<<(poly lhs, U t) noexcept {
		return lhs <<= t;
	}
	constexpr poly Rev() const noexcept {
		return poly{vector<T>::rbegin(), vector<T>::rend()};
	}
	constexpr poly &shrink() {
		if (!vector<T>::empty())
			for (u32 k = deg() - 1; k; --k)
				if (f[k] != 0) return redeg(k + 1);
		return *this;
	}

   private:
	constexpr poly &ClearL(u32 degree) noexcept {
		return std::fill_n(vector<T>::begin(), degree >> 1, 0), *this;
	}
	constexpr poly &ClearH(u32 degree) noexcept {
		return std::fill_n(vector<T>::begin() + (degree >> 1), degree >> 1, 0),
			   *this;
	}
	constexpr poly &InvHelp(poly p, poly q, u32 degree) noexcept {
		times(p, q, degree).ClearL(degree);
		times(p, q, degree), redeg(degree);
		for (u32 k = degree >> 1; k < degree; ++k) f[k] = -p[k];
		return *this;
	}

   public:
	constexpr poly Inv(u32 degeree) const noexcept {
		poly ret{1 / f.front()};
		for (u32 l = 2; l < degeree * 2; l *= 2)
			ret.InvHelp(slice(l), ret.slice(degeree), l);
		return ret.redeg(degeree);
	}
	constexpr poly Div(poly g, u32 degree) const noexcept {
		if (deg() == 0) return {};
		u32 t = toLim(degree);
		poly f0 = slice(t >> 1), g0 = g.Inv(t >> 1);
		poly q = times(f0, g0, t).slice(t >> 1);
		times(q, g, t).ClearL(t);
		for (u32 k = t >> 1; k < Utility::min(t, deg()); ++k) q[k] -= f[k];
		times(q, g0, t), f0.redeg(degree);
		for (u32 k = t >> 1; k < Utility::min(t, degree); ++k) f0[k] = -q[k];
		return f0;
	}
	constexpr poly Sqrt(u32 degree) const noexcept {
		poly ret{1};
		for (u32 l1 = 2, l2 = 4; l1 < degree * 2; l1 = l2, l2 <<= 1)
			(ret = (ret + slice(l1).Div(ret, l1)) * (T(2).inv())).redeg(l1);
		return ret.redeg(degree);
	}
	constexpr poly Ln(u32 degree) const noexcept {
		return Dvt(degree).Div(slice(degree), degree).Igt(degree);
	}
	constexpr poly Exp(u32 degree) const noexcept {
		poly ret{1};
		for (u32 l1 = 2, l2 = 4; l1 < degree * 2; l1 = l2, l2 <<= 1)
			ret *= poly{1} - ret.Ln(l1) + slice(l1);
		return ret.redeg(degree);
	}
	constexpr poly Pow(u32 degree, T pw) const noexcept {
		return (Ln(degree) * pw).Exp(degree);
	}
	friend constexpr poly operator/(poly p, const poly &q) noexcept {
		if (p.deg() < q.deg()) p.redeg(q.deg());
		return p.Rev().Div(q.Rev(), p.deg() - q.deg() + 1).Rev();
	}
	friend constexpr poly operator%(const poly &p, const poly &q) noexcept {
		poly g = p / q, h = p - g * q;
		return h.shrink();
	}
	constexpr std::pair<poly, poly> Divt(const poly &g) const noexcept {
		poly p = f / g, q = f - g * p;
		return {p, q.shrink()};
	}
	constexpr poly ChirpZ(T c, u32 m) const noexcept {
		poly p(deg(), [&](u32 k) { return f[k] / c.pow(k * (k - 1ll) >> 1); }),
			q(deg() + m - 1, [&](u32 k) { return c.pow(k * (k - 1ll) >> 1); });
		std::reverse(p.begin(), p.end()), p *= q;
		for (u32 k = 0; k < m; ++k)
			p[k] = p[k + deg() - 1] / c.pow((k * (k - 1ll)) >> 1);
		return p.redeg(m);
	}
#undef f
};
template <typename T, typename _Trans>
_Trans poly<T, _Trans>::trans;

enum FWT_type { Or = 0, And = 1, Xor = 2 };
template <FWT_type type>
struct FWT {};
template <>
struct FWT<Or> {
	template <Trans_mode mode, typename Iterator>
	constexpr void trans(Iterator f, u32 n) const noexcept {
		for (u32 k = 1; k < n; k <<= 1)
			for (u32 i = 0; i < n; ++i)
				if (i & k) {
					if constexpr (mode == DIF)
						f[i] += f[i ^ k];
					else
						f[i] -= f[i ^ k];
				}
	}
};
template <>
struct FWT<And> {
	template <Trans_mode mode, typename Iterator>
	constexpr void trans(Iterator f, u32 n) const noexcept {
		for (u32 k = 1; k < n; k <<= 1)
			for (u32 i = 0; i < n; ++i)
				if (i & k) {
					if constexpr (mode == DIF)
						f[i ^ k] += f[i];
					else
						f[i ^ k] -= f[i];
				}
	}
};
template <>
struct FWT<Xor> {
	template <typename T>
	static constexpr T Inv2 = T(1) / 2;
	template <Trans_mode mode, typename Iterator>
	static constexpr void trans(Iterator f, u32 n) noexcept {
		auto inv2 = Inv2<std::remove_reference_t<decltype(f[0])>>;
		for (u32 p = 1, l = 2; p < n; p = l, l <<= 1)
			for (u32 t = 0; t < n; t += l)
				for (u32 k = 0; k < p; ++k) {
					auto x = f[t | k], y = f[t | p | k];
					if constexpr (mode == DIF)
						f[t | k] += y, f[t | p | k] = x - y;
					else
						(f[t | k] += y) *= inv2, f[t | p | k] = (x - y) * inv2;
				}
	}
};
template <typename Trans, typename Iterator>
constexpr void times(Iterator f, Iterator g, u32 deg, Trans trans = Trans{}) {
	trans.template trans<DIF>(f, deg);
	trans.template trans<DIF>(g, deg);
	dot(f, g, deg);
	trans.template trans<DIT>(f, deg);
	return f;
}

}  // namespace Poly

#undef UseStd
#undef OPERATOR
}  // namespace Default_Source
using namespace Default_Source;

constexpr i64 N = 1e15;
IO::InStream<0> qin;
using Arithmetic::Z64;
Z64<0> f[14], g[14];
struct pair {
	i64 x;
	Z64<0> val;
	pair() = default;
	pair(i64 _x, Z64<0> _val) : x{_x}, val{_val} {}
	bool operator<(const pair &ot) const noexcept { return x < ot.x; }
};
std::vector<i64> pre;
std::vector<pair> pref, preg;
Z64<0> calc(i64 n, const std::vector<pair> &pre_func) {
	return (--std::upper_bound(pre_func.begin(), pre_func.end(), pair(n, {})))
		->val;
}
Z64<0> get(i64 n, i32 t, const std::vector<pair> &pre_func) {
	return n % t ? Z64<0>(0) : calc(n / t, pre_func);
}

signed main() {
	auto [f0, g0, T, mod] = qin.reads<i32, i32, i32, i64>();
	Z64<0>::set(mod), f[0] = f0, g[0] = g0;
	for (i32 n = 1; n <= 13; ++n)
		f[n] = max<i64>(Z64<0>(n).val(),
						(g[n / 2] + g[n / 3] + g[n / 5] + g[n / 7]).val()),
		g[n] = max<i64>(Z64<0>(n).val(),
						(f[n / 2] + f[n / 3] + f[n / 4] + f[n / 5]).val());
	pre.reserve(84078), pre.emplace_back(0);
	for (i32 n : {1, 2, 3, 5, 7, 9, 11, 13})
		for (i64 _2 = 1, __2 = N / n; _2 <= __2; _2 <<= 1)
			for (i64 _3 = 1, __3 = N / (n * _2); _3 <= __3; _3 *= 3)
				for (i64 _5 = 1, __5 = N / (n * _2 * _3); _5 <= __5; _5 *= 5)
					for (i64 _7 = 1, __7 = N / (n * _2 * _3 * _5); _7 <= __7;
						 _7 *= 7)
						pre.emplace_back(n * _2 * _3 * _5 * _7);
	std::sort(pre.begin(), pre.end()),
		pre.erase(std::unique(pre.begin(), pre.end()), pre.end());
	for (i32 n = 13; n >= 1; --n) f[n] -= f[n - 1], g[n] -= g[n - 1];
	pref.reserve(pre.size() + 1), preg.reserve(pre.size() + 1);
	for (auto n : pre)
		if (n <= 13)
			pref.emplace_back(n, f[n]), preg.emplace_back(n, g[n]);
		else {
			Z64<0> posf = get(n, 2, preg) + get(n, 3, preg) + get(n, 5, preg) +
						  get(n, 7, preg),
				   posg = get(n, 2, pref) + get(n, 3, pref) + get(n, 4, pref) +
						  get(n, 5, pref);
			pref.emplace_back(n, posf), preg.emplace_back(n, posg);
		}
	for (u32 t = 1; t < pre.size(); ++t)
		pref[t].val += pref[t - 1].val, preg[t].val += preg[t - 1].val;
	for (i64 n; T; --T)
		qin.read(n),
			printf("%lld %lld\n", calc(n, pref).val(), calc(n, preg).val());

	return 0;
}

Details

answer.code: In function ‘constexpr T Default_Source::Arithmetic::mul(T, T, T)’:
answer.code:272:47: error: ‘u128’ was not declared in this scope
  272 |         else if constexpr (sizeof(T) < sizeof(u128))
      |                                               ^~~~