https://agc013.contest.atcoder.jp/tasks/agc013_e
一応知ってるテクニックの組み合わせにはなるんだけど…難しいですね。
Degwerさんの数え上げテクニック集で紹介されている問題です。
とりあえず自明DPを考えると
dp(x)=Σ(y=1...x)y^2dp(x-y)となります。こういうy=1からのΣの式はy=1の場合を分離するとうまく行くことが多いです。蟻本にも載ってますね。
そのように式変形すれば、
dp(x)=dp(x-1)+Σ(y=1...x-1)(y+1)^2dp(x-y)=dp(x-1)+Σ(y=1...x-1)y^2dp(x-y)+2Σ(y=1...x-1)ydp(x-y)+Σ(y=1...x-1)dp(x-y)
=2*dp(x-1)+2dp'(x-1)+dp''(x-1)(ただしdp'(x-1)=Σydp(x-y), dp''(x-1)=Σdp(x-y))
となって、dp'(x),dp''(x)もdp(x-1),dp'(x-1),dp''(x-1)の線型結合で表せるので行列累乗が使えます。
おいてはいけない点についても掛ける行列を多少変えれば対応できてO(MlogN)です。
Editorialでは場合の数に帰着していますが思いつくわけがありませんね。でもこういう和を場合の数に帰着して計算することがあるのは勉強になりました。
int N, M, C; int A[MAX_N]; int L[MAX_N], R[MAX_N]; void add_interval(int l, int r) { if(l <= r) { L[C] = l; R[C] = r; C++; } } mat skip(ll k) { mat res(3, vl(3, 0)); res[0][0] = (k + 1) * (k + 1) + 1; res[0][1] = 2 * (k + 1); res[0][2] = (k + 1) * (k + 1); res[1][0] = k + 1; res[1][1] = 1; res[1][2] = k + 1; res[2][0] = 1; res[2][1] = 0; res[2][2] = 1; rep(i, 0, 3) rep(j, 0, 3) res[i][j] %= mod; return res; } void solve() { cin >> N >> M; rep(i, 0, M) cin >> A[i]; C = 0; mat BB = skip(0); if(M != 0) { add_interval(0, A[0] - 1); rep(i, 0, M - 1) { add_interval(A[i] + 1, A[i + 1] - 1); } add_interval(A[M - 1] + 1, N); mat D(3, vl(3, 0)); D[0][0] = 1; D[1][1] = 1; D[2][2] = 1; rep(i, 0, C) { D = mat_mul(mat_pow(BB, R[i] - L[i]), D); if(i != C - 1) { D = mat_mul(skip(L[i + 1] - R[i] - 1), D); } } cout << D[0][2] << "\n"; } else { cout << mat_pow(BB, N)[0][2] << "\n"; } }