提交记录 17515


用户 题目 状态 得分 用时 内存 语言 代码长度
Misaka 1004. 【模板题】高精度乘法 Accepted 100 35.934 ms 28332 KB C++ 6.61 KB
提交时间 评测时间
2022-03-21 12:41:02 2022-03-21 12:41:07
#if !IN_LOCAL
#define ReportToStdErr(...) __VA_ARGS__
#endif // !IN_LOCAL
#include<cstdio>
#include<cmath>
#include<cstring>
constexpr const double Pi = 3.141592653589793115997963468544185161590576171875;
struct Complex
{
	double real, imag;
	constexpr Complex(const double _real = 0, const double _imag = 0) :real(_real), imag(_imag) {}
	inline constexpr Complex operator+(const Complex& other)const
	{
		return Complex(real + other.real, imag + other.imag);
	}
	inline constexpr Complex& operator+=(const Complex& other)
	{
		real += other.real;
		imag += other.imag;
		return *this;
	}
	inline constexpr Complex operator-()const
	{
		return Complex(-real, -imag);
	}
	inline constexpr Complex operator-(const Complex& other)const
	{
		return Complex(real - other.real, imag - other.imag);
	}
	inline constexpr Complex operator*(const Complex& other)const
	{
		return Complex(real * other.real - imag * other.imag, real * other.imag + imag * other.real);
	}
	inline constexpr Complex& operator*=(const Complex& other)
	{
		return *this = *this * other;
	}
};
constexpr const size_t sizeof_Complex = sizeof(Complex);
inline constexpr Complex Conjugate(const Complex& val)
{
	return Complex(val.real, -val.imag);
}
inline constexpr Complex AddOne(const Complex& val)
{
	return Complex(val.real + 1, val.imag);
}
inline constexpr Complex DivFour(const Complex& val)
{
	return Complex(val.real / 4, val.imag / 4);
}
char buf_in[2097152], buf_out[2097152], * p_buf_out = buf_out;
const char* p_buf_in = buf_in;
constexpr const char& buf_in_0 = buf_in[0], & buf_in_1 = buf_in[1], & buf_in_2 = buf_in[2], * const buf_inP1 = buf_in + 1, * const buf_inP2 = buf_in + 2, * const buf_inP3 = buf_in + 3,
* const buf_inP4 = buf_in + 4;
int a[524288];
Complex f[524288], g[524288], h[524288], unit_root[524288], weight[524288];
constexpr int& a_0 = a[0], & a_1 = a[1];
constexpr Complex& f_0 = f[0], & g_0 = g[0], & h_0 = h[0], & unit_root_1 = unit_root[1], * const hP1 = h + 1;
int main()
{
	fread(buf_in, 1, 2097152, stdin);
	for (; *p_buf_in >= '0'; ++p_buf_in);
	int len_a = 0;
	const char* p = p_buf_in;
	for (; p >= buf_inP4; a[len_a++] = (*(p - 1) ^ '0') + (*(p - 2) ^ '0') * 10 + (*(p - 3) ^ '0') * 100 + (*(p - 4) ^ '0') * 1000, p -= 4);
	if (p == buf_inP3)
	{
		a[len_a++] = (buf_in_0 ^ '0') * 100 + (buf_in_1 ^ '0') * 10 + (buf_in_2 ^ '0');
	}
	else if (p == buf_inP2)
	{
		a[len_a++] = (buf_in_0 ^ '0') * 10 + (buf_in_1 ^ '0');
	}
	else if (p == buf_inP1)
	{
		a[len_a++] = buf_in_0 ^ '0';
	}
	for (int i = 0; i < len_a; f[i >> 1].real = a[i], f[i >> 1].imag = a[i | 1], i += 2);
	int length = (len_a >> 1) + (len_a & 1) - 1, len_ans = len_a - 1;
	for (; *p_buf_in < '0'; ++p_buf_in);
	const char* const start_pos = p_buf_in, * const start_posP4 = start_pos + 4;
	for (; *p_buf_in >= '0'; ++p_buf_in);
	for (len_a = 0; p_buf_in >= start_posP4; a[len_a++] = (*(p_buf_in - 1) ^ '0') + (*(p_buf_in - 2) ^ '0') * 10 + (*(p_buf_in - 3) ^ '0') * 100 + (*(p_buf_in - 4) ^ '0') * 1000, p_buf_in -= 4);
	if (p_buf_in == start_pos + 3)
	{
		a[len_a++] = (*start_pos ^ '0') * 100 + (*(start_pos + 1) ^ '0') * 10 + (*(start_pos + 2) ^ '0');
	}
	else if (p_buf_in == start_pos + 2)
	{
		a[len_a++] = (*start_pos ^ '0') * 10 + (*(start_pos + 1) ^ '0');
	}
	else if (p_buf_in == start_pos + 1)
	{
		a[len_a++] = *start_pos ^ '0';
	}
	for (int i = a[len_a] = 0; i < len_a; g[i >> 1].real = a[i], g[i >> 1].imag = a[i | 1], i += 2);
	const int lengthD2 = (length = ((length += (len_a >> 1) + (len_a & 1)) > 1 ? 1 << ((length = 31 - __builtin_clz(length) + ((length & (length - 1)) > 0))) : 2)) >> 1;
	len_ans += len_a;
	unit_root_1 = 1;
	double angle = Pi;
	for (int i = 2; i < length;)
	{
		angle /= 2;
		const Complex gn(cos(angle), sin(angle));
		for (int j = i, end_j = i *= 2; j < end_j; unit_root[j | 1] = (unit_root[j] = unit_root[j >> 1]) * gn, j += 2);
	}
	for (int i = lengthD2, end_j = 1; i; i >>= 1, end_j *= 2)
	{
		const Complex* const start_pur = unit_root + i;
		Complex* p1 = f, * p2 = g;
		for (int j = 0; j < end_j; ++j, p1 += i, p2 += i)
		{
			const Complex* p_unit_root = start_pur;
			for (int k = 0; k < i; ++k, ++p1, ++p2, ++p_unit_root)
			{
				const Complex t1 = *p1, t2 = *p2;
				Complex& t3 = *(p1 + i), & t4 = *(p2 + i);
				*p1 += t3;
				*p2 += t4;
				t3 = (t1 - t3) * *p_unit_root;
				t4 = (t2 - t4) * *p_unit_root;
			}
		}
	}
	(h_0 = f_0 * g_0).real += 2 * f_0.imag * g_0.imag;
	memcpy(weight, unit_root + lengthD2, lengthD2 * sizeof_Complex);
	for (int i = 0; i < lengthD2; weight[i + lengthD2] = -weight[i], ++i);
	for (int i = 0, p = 0; i < length; ++i)
	{
		if (p < i)
		{
			const Complex t = weight[p];
			weight[p] = weight[i];
			weight[i] = t;
		}
		for (int j = lengthD2; (p ^= j) < j; j >>= 1);
	}
	for (int i = 1; i < lengthD2;)
	{
		const int end_j = i * 2, t = i + end_j - 1;
		for (int j = i; j < end_j; h[j] = f[j] * g[j] - DivFour(AddOne(weight[j]) * (f[j] - Conjugate(f[t - j])) * (g[j] - Conjugate(g[t - j]))), ++j);
		i = end_j;
	}
	for (int i = lengthD2, t = lengthD2 + length - 1; i < length; h[i] = f[i] * g[i] - DivFour(AddOne(weight[i]) * (f[i] - Conjugate(f[t - i])) * (g[i] - Conjugate(g[t - i]))), ++i);
	for (int i = 1, end_j = lengthD2; i < length; i *= 2, end_j >>= 1)
	{
		const Complex* const start_pur = unit_root + i;
		Complex* p = h;
		for (int j = 0; j < end_j; ++j, p += i)
		{
			const Complex* p_unit_root = start_pur;
			for (int k = 0; k < i; ++k, ++p, ++p_unit_root)
			{
				const Complex t = *(p + i) * *p_unit_root;
				*(p + i) = *p - t;
				*p += t;
			}
		}
	}
	const double inv_length = 1.0 / length;
	unsigned long long up = static_cast<long long>(h_0.real * inv_length + 0.5);
	a_0 = up % 10000;
	up /= 10000;
	a_1 = (up += static_cast<long long>(h_0.imag * inv_length + 0.5)) % 10000;
	up /= 10000;
	for (int i = 2, p = 0; i <= len_ans; i += 2)
	{
		a[i] = (up += static_cast<long long>(h[length - (i >> 1)].real * inv_length + 0.5)) % 10000;
		up /= 10000;
		a[i | 1] = (up += static_cast<long long>(h[length - (i >> 1)].imag * inv_length + 0.5)) % 10000;
		up /= 10000;
	}
	for (; len_ans && a[len_ans] == 0; --len_ans);
	unsigned short t = a[len_ans];
	if (t > 9)
	{
		if (t > 99)
		{
			if (t > 999)
			{
				*(p_buf_out++) = static_cast<char>(t / 1000) ^ '0';
				t %= 1000;
			}
			*(p_buf_out++) = static_cast<char>(t / 100) ^ '0';
			t %= 100;
		}
		*(p_buf_out++) = static_cast<char>(t / 10) ^ '0';
		t %= 10;
	}
	for (*(p_buf_out++) = static_cast<char>(t) ^ '0'; len_ans--;)
	{
		*(p_buf_out++) = static_cast<char>((t = a[len_ans]) / 1000) ^ '0';
		*(p_buf_out++) = static_cast<char>((t %= 1000) / 100) ^ '0';
		*(p_buf_out++) = static_cast<char>((t %= 100) / 10) ^ '0';
		*(p_buf_out++) = static_cast<char>(t % 10) ^ '0';
	}
	fwrite(buf_out, 1, p_buf_out - buf_out, stdout);
	return 0;
}

CompilationN/AN/ACompile OKScore: N/A

Testcase #135.934 ms27 MB + 684 KBAcceptedScore: 100


Judge Duck Online | 评测鸭在线
Server Time: 2024-12-05 10:17:17 | Loaded in 0 ms | Server Status
个人娱乐项目,仅供学习交流使用 | 捐赠