본문 바로가기

분할 정복을 이용한 거듭제곱/행렬 거듭제곱

백준 31987번: ESC와 쿼리

https://www.acmicpc.net/problem/31987


문제에서 요구하는 것은 다음과 같습니다.

$$ \sum_{x=i}^{j} a_{kx} + b_{kx} + c_{kx} \Longleftarrow \frac{d^{n}}{dx^{n}} \bigg ( e^{x}sin(x)cos(x) \bigg ) = a_{n} e^{x}sin^{2}(x) + b_{n} e^{x} cos^{2} (x) + c_{n} e^{x} sin(x) cos(x) $$

 

그런데 $e^{x}$ 는 미분해도 자기 자신이며, $sin(x)$ 를 미분하면 $cos(x)$ 가 되고, $cos(x)$ 를 미분하면 $-sin(x)$ 가 된다는 것이 알려져 있습니다. 따라서, 이 값을 구하기 위해 초항부터 구해나가며 관찰을 해도 좋지만, $n$ 번 미분된 식을 한 번 더 미분해도 되겠습니다.

 

$$ \frac{d^{n+1}}{dx^{n+1}} \bigg ( a_{n}e^{x}sin^{2}(x) + b_{n}e^{x}cos^{2}(x) + c_{n}e^{x}sin(x)cos(x) \bigg ) $$

$$ = \big ( a_{n} - c_{n} \big ) e^{x}sin^{2}(x) + \big ( b_{n}+c_{n}  \big ) e^{x}cos^{2}(x) + \big ( 2a_{n} - 2b_{n} + c_{n} \big ) e^{x}sin(x)cos(x) $$

$$ \equiv a_{n+1}e^{x}sin^{2}(x) + b_{n+1}e^{x}cos^{2}(x) + c_{n+1}e^{x}sin(x)cos(x) $$

 

마지막 줄로부터, 점화식을 행렬로 쉽게 표현할 수 있습니다.

$$ \begin{bmatrix} a_{n} \\ b_{n} \\ c_{n} \end{bmatrix} = \begin{bmatrix} 1 & 0 & -1 \\ 0 & 1 & 1 \\ 2 & -2 & 1 \end{bmatrix}^{n-1} \begin{bmatrix} a_{1} \\ b_{1} \\ c_{1} \end{bmatrix} $$

 

그런데, 구해야 하는 식은 $\sum_{x=i}^{j} \big ( a_{kx} + b_{kx} + c_{kx} \big ) $ 입니다. 나이브하게 구하게 되면 쿼리 당 $O(3^{3} \times \log ik + 3^{3} \times (j-i))$ 라는, 너무 큰 복잡도가 나오기 때문에, 식을 간소화하고 복잡도를 낮출 방법을 찾아보겠습니다.


제일 먼저, 위 점화식으로부터

$$ a_{n+1} + b_{n+1} = a_{n} + b_{n} $$

임을 알 수 있으며, 초항을 구해보면 이는 0입니다! 그러므로 구해야 하는 식은

$$ \sum_{x=i}^{j} c_{kx} $$

로 간소화됩니다. 추가적으로, $a_{n} + b_{n} = 0$ 으로부터 점화식도 간소화됩니다.

$$ \begin{bmatrix} a_{n} \\ c_{n} \end{bmatrix} = \begin{bmatrix} 1 & -1 \\ 4 & 1 \end{bmatrix}^{n-1} \begin{bmatrix} -1 \\ 1 \end{bmatrix} $$

하지만, 여전히 복잡도의 차수는 동일하기 때문에, 거듭제곱의 합을 빠르게 구할 방법을 찾아야만 합니다. 행렬의 등비수열로 접근한다면 예제까지는 충분히 답을 얻을 수 있으나, 역행렬이 보장되게 만드는 게 쉽진 않을 것 같습니다. 그렇다면, 일반항 밖에 선택의 여지가 없습니다.


Wolframalpha, 또는 직접 풀어 점화식의 일반항을 얻으면 아마 다음과 같은 형태로 나타날 것입니다. 직접 푸는 과정은 나중에 첨부하도록 하겠습니다.

$$ c_{n} = \frac{(1+2i)^{n} + (1-2i)^{n}}{2} $$

이 식을 다양하게 변형할 수 있으나, 문제의 답을 구할 때는 등비수열의 합으로 구하는 게 가장 자연스러울 것 같습니다.

$$ \sum_{x=i}^{j} c_{kx} = \frac{1}{2} \times \sum_{x=i}^{j} \bigg ( (1+2i)^{kx} + (1-2i)^{kx} \bigg ) $$

$$ = \frac{1}{2} \times \left [ \left \{ \sum_{x=0}^{j} (1+2i)^{kx} - \sum_{x=0}^{i-1} (1+2i)^{kx} \right \} + \left \{ \sum_{x=0}^{j} (1-2i)^{kx} - \sum_{x=0}^{i-1} (1-2i)^{kx} \right \} \right ] $$

$$ = \frac{1}{2} \times \left \{ \frac{(1+2i)^{k(j+1)} - (1+2i)^{ki}}{(1+2i)^{k} - 1} + \frac{(1-2i)^{k(j+1)} - (1-2i)^{ki}}{(1-2i)^{k} - 1} \right \} $$

놀랍게도, 이 복잡한 식에서 통분을 하게 되면 상당히 예쁜 식이 나오게 됩니다.

$$ \sum_{x=i}^{j} c_{kx} = \frac{5^{k} \left ( c_{kj} - c_{k(i-1)} \right ) - c_{k(j+1)} + c_{ki}}{5^{k} + 1 - 2 \times c_{k}} $$

이 식을 바탕으로 문제를 풀게 되면 시간 복잡도는 $O(Q \times 2^{3} \times \log jk)$ 입니다.

 

C/C++은 시간제한이 0.2s이기 때문에, 이 정도의 복잡도를 가지고도 구현에서 최적화를 하지 않으면 시간 초과될 여지가 있습니다. 제 코드는 다음과 같습니다.

#include <iostream>
#define mod 1000000007

typedef unsigned long long ull;

ull power_modular (ull X, ull Y) {
	ull result = 1;
	while(Y) {
		if (Y % 2) {
			result = result * X % mod;
		}
		X = X * X % mod;
		Y /= 2;
	}
	return result;
}

typedef struct {
	ull array[2][2];
} Matrix;

Matrix matrix_multiply_modular (Matrix A, Matrix B) {
	Matrix result = {{{}}};
	result.array[0][0] = (A.array[0][0] * B.array[0][0] + A.array[0][1] * B.array[1][0]) % mod;
	result.array[0][1] = (A.array[0][0] * B.array[0][1] + A.array[0][1] * B.array[1][1]) % mod;
	result.array[1][0] = (A.array[1][0] * B.array[0][0] + A.array[1][1] * B.array[1][0]) % mod;
	result.array[1][1] = (A.array[1][0] * B.array[0][1] + A.array[1][1] * B.array[1][1]) % mod;
	return result;
}

Matrix matrix_power_modular (Matrix A, long long P) {
	Matrix result = {{{1,0},{0,1}}};
	while (P) {
		if (P % 2) {
			result = matrix_multiply_modular(result, A);
		}
		A = matrix_multiply_modular(A, A);
		P /= 2;
	}
	return result;
}

using namespace std;
int main() {
	ios_base::sync_with_stdio(false);
	cin.tie(nullptr);
	cout.tie(nullptr);
	
	int Q;
	cin >> Q;
	while(Q--)
	{
		long long i, j, k;
		cin >> i >> j >> k;

		Matrix base = {{{1, mod - 1},{4, 1}}};
		
		Matrix temp1 = matrix_power_modular(base, k);
		Matrix temp2 = matrix_power_modular(base, k * (i - 1));
		Matrix temp3 = matrix_multiply_modular(temp1, temp2);
		Matrix temp4 = matrix_power_modular(base, k * j);
		Matrix temp5 = matrix_multiply_modular(temp1, temp4);
		
		ull temp = power_modular(5, k);
		
		ull divisor	= (temp + 1 + (mod - 2) * temp1.array[1][1]) % mod;
		ull dividend = (temp * ((temp4.array[1][1] + (mod - 1) * temp2.array[1][1]) % mod) + (mod - 1) * temp5.array[1][1] + temp3.array[1][1]) % mod;
		
		ull answer = dividend * power_modular(divisor, mod - 2) % mod;
		cout << answer << "\n";
	}

	return 0;
}

(C++20, 112ms, 2020KB, 제출번호 86595282)